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Introduction 


In Unit 9, Section 5, we discussed the solution of sets of homogeneous 
simultaneous equations. In this unit, we examine a related problem which often 
arises in practice. For instance, in Unit 24, you will be studying systems of 
particles and springs like those shown in Figure 1. 


Y 


afm Fm 
ri 
# 


Figure 1 


The equations of motion for this system are a set of linear simultaneous 
differential equations. In Unit 22, you will see that a first step in solving such a set 
of differential equations is to find the column vectors x which satisfy a set of 
algebraic linear equations of the form 


Ax = (x, (1) 


where A is a square matrix which is determined by the set of differential 

equations and 4 is a number whose value is to be determined. So, before we can 
do the work in Units 22 and 24, we first have to discuss methods of solution of 
Equation (1). In general the only solution is x = 0, but there are particular values 
of 4 for which solutions other than x = 0 are possible. These particular values of A, 
and the associated non-zero column vectors x, enable us to solve the equations of 
motion for the system shown in Figure 1. 


Equations of the form Ax = 4x also arise in a large variety of other mechanical 
and electrical systems in engineering. There are also statistical, numerical and 
other non-engineering situations in which such equations arise. In fact, they arise 
so frequently that it is worth discussing the problem in its own right, which is 
what we will do in this unit. The values of 2 for which Equation (1) has non-zero 
solutions are called the eigenvalues of A, and the non-zero solutions are called the 
eigenvectors of A. 


Study guide 


You should study Section 1 and Subsections 2.1 and 2.2 before watching the 
television programme. 


Section 4 describes a computer package which you can use to find eigenvalues and 
eigenvectors of a given matrix. You might like to use this package when studying 
Units 22 and 24, as well as in this unit. The exercises at the end of Section 4 are 
there for use at the computer terminal. You are expected to do as many as you 
have time for. However, you may not have time to do them all, so start with the 
ones that interest you the most. 


The exercises at the end of Sections 1, 2 and 3 are there to give you extra practice. 
It is not necessarily expected that you will do them when first reading this unit. 
You might like to use them for revision purposes. However, you are encouraged to 
do the end of unit test in Section 6 immediately after studying the unit. 


MST204 21.1 5 


| The theoretical eigenvalue problem 


1.1 Introduction 


In this section, we discuss a method for finding non-zero vectors x which satisfy 
the matrix equation 


Ax = 4x. (1) 


Here, A is a known square matrix, but at the beginning of the problem, 4 is an 
unknown number. 


In fact, this matrix equation turns out to be a set of homogeneous linear 
simultaneous equations, as the following example demonstrates. 


Example 1 


3 2 
A= 5 A! and x = ay show that the matrix equation Ax = Ax may be 
S 


written as a set of homogeneous linear simultaneous equations. 


Solution 
The matrix equation is 


2 sil-n} 
2 3 |] x2 X> 
This gives us the equations 

Jhy oe 2%5-= AX, 

Pie he = AX, 


and these can be re-written as a pair of homogeneous linear simultaneous 
equations 


(5 — A)x, + Lis GQ 
2x, + (5 — A)x, = 0. 


More generally, the equation 


Ax = Ax 
can be re-written as 
Ax = dIx, 
SO 
Ax — AIx = 0, 
or 
(A — AI)x = 0. (2) 


The matrix equation (2) is just the matrix form of a set of homogeneous linear 
simultaneous equations. Such equations were discussed in Section 5 of Unit 9. You 
will remember that homogeneous sets of equations always have the ‘trivial 
solutions’ x = 0. However, if the equations are linearly dependent they will have an 
infinite number of other solutions as well. Now, linear dependence of 
homogeneous sets of equations is guaranteed if the determinant of the left-hand- 
side coefficients is zero. So, there are non-trivial solutions to the equations 


Ax = Ax 
ie (hk — de 0 (2) 
if A satisfies 


det (A — AI) = 0. (3) 


I is the identity matrix defined 
in Unit 20, Subsection 4.1. 


For this application of 
determinants see Unit 20, 
Subsection 5.3. 
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Example Z 
The equations from Example 1: 
(5 — A)x, + 2X, = 0 
2x, + (5 — A)x, = 0 
have non-trivial solutions if 
5—A 2 


2 Ae ees 


From the determinant in this example, we could have gone on to find /. 
In general, the determinant condition (3) can be used to find values of A for which 
the equations have nontrivial solutions. These values of 4 are called eigenvalues of A. 


For each eigenvalue, we get a particular set of equations, which can be solved for 
x. These solutions are called eigenvectors, and the whole problem is often referred 
to loosely as the eigenvalue problem. 


The rest of this section is made up mostly of examples and exercises. This is to 
help you gain some experience in solving simple eigenvalue problems. 


1.2 The 2 x 2 eigenvalue problem 


Example 3 
Find the eigenvalues of the matrix 


Solution 


The obvious way to approach this is to start with the definition of an eigenvalue, 
and to find the values of 2 for which the equation Ax = Ax, that is 


5 2]) X1] _ 4 Pl 
2 5 X72 = X2 : 
has non-trivial solutions. 


We know from Example 1 that this matrix equation gives us the separate 
equations 
(S = A)x4 5 2X> se () 
2x, + (5 — A)x, = 0. 


For non-trivial solutions, we require the determinant of the left-hand-side 
coefficients to be zero. So, 


5—A Z 
2 oes 
This gives us the equation 
(5 — A)? —27=0 (4) 
SO ((5 — A) — 2)((5 — A) + 2) =0 
or (3 — A)(7 — A) = 0. 


Thus we can find two values, 2 = 3 or 4 = 7, for which the equations have non- 
trivial solutions. Hence, the eigenvalues of A are 3 and 7. 


Each time we use the above method to find the eigenvalues of a matrix A, we have 
to find the values of 4 which satisfy det(A — AI) = 0. So, it is normal practice to 
short-cut the work shown by going directly from the matrix A to the determinant 
equation det (A — AI) = 0. For instance, in the example above we could easily have 
written down 


5-A Zz 
2 5—A 
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straight away, and solved the resulting equation 
(5—AyY —27=0 (4) 
as before. 


The equation which we obtain by expanding the determinant (such as Equation 
(4) above) is often called the characteristic equation of the matrix. 


Example 4 
Find the eigenvalues of 


Solution 
The condition det (A — AI) = 0 gives 
2-—A 1 


= 0 
1 2—A 


and hence we obtain the characteristic equation 
(2—A)y? —17 =0 

Le. (3 —A)(1 — A) = 0. 

Hence, the eigenvalues of A are 3 and 1. 


In general, the method for finding the eigenvalues of a 2 x 2 matrix is as follows. 


Procedure 1.2(a): To find the eigenvalues of a 2 x 2 matrix 


To find the eigenvalues of 


write down the characteristic equation det (A — AI) = 0, that is, 


a4, —A 
a21 


and solve for 2. 


In theory, this method can be extended to find the eigenvalues of any square 
matrix A. 


Exercise 1 
Find the eigenvalues of the matrices 


(i) | 7 3 (ii) | 3 2 
3 7 1 4 
[Solution on p. 45 | 


Once the eigenvalues of a matrix have been found we can go on to find the 
corresponding eigenvectors. 


Example 5 
Find the eigenvectors of the matrix 


Solution 


We know from Example 3 that the eigenvalues are 3 and 7. To find the 
corresponding eigenvectors we consider each eigenvalue in turn and solve 
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(A — Al)x =0 (2) 
for x. 
Case 1 = 3: in this case Equation (2) becomes 
(5 — 3)x,; + 2x, =Q 
2x, + (5 — 3)x, = 0. 
The fact that both these equations give us the same information 
x, +x, =0 


is not surprising, as the determinant condition forces the equations to be linearly 
dependent. 


To solve these equations we put x, = k, where k is an arbitrary number. 
Substituting this into the equation above, we get x, = —k. In this way we obtain 
an infinite number of solutions of the form 


So, with 1 = 3, any vector of the form 


—1 
= i 
a 
will satisfy Equation (2). 


When we give an eigenvector corresponding to a given eigenvalue, it is common 
practice to omit the k, and just to say that an eigenvector corresponding to the 


eigenvalue 3 is a In doing this, it is understood that B ih a etc. would 
be equally correct, and that an eigenvector is not unique. 
Case 1 = 7: in this case Equation (2) becomes 
(5 — 7)x, + 2x, =0 
2x, + (5 — 7)x, = 0. 
These give us the single equation 


— X41 + X2 == 0). 


, k 
Putting x, = k into this equation, we get x, = k. Thus, any vector of the form i 


is a solution of Ax = 7x and so an eigenvector corresponding to the eigenvalue 7 
ee 
is 


Procedure 1.2(b): To find the eigenvectors of a 2 x 2 matrix 


To find the eigenvectors of 


1. Find the eigenvalues 4, and A, of A. 


2. For each of these values of 4 solve (A — AI)x = 0 for x. That is, 
solve 


(a4; — Ax)X1 + Ay2X2 =0 
Az1X1 + (Ap2 — A,)X2 = 9 


{or R= 1,2. 


Then any non-trivial solution to these equations will give an 
eigenvector corresponding to the eigenvalue A,. 


Again, the method can in theory be extended to find the eigenvectors of any square 
matrix A. 


This approach is first 
discussed in Section 2 of 
Unit 9. 
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Exercise 2 
Use the eigenvalues you found in Exercise 1 to find eigenvectors of 


(i) k a. my is 
3 , eS 


[Solution on p. 45] 


1.3. The 3 x 3 eigenvalue problem 


The problem of finding eigenvalues and eigenvectors of a 3 x 3 matrix can be 
tackled in the same way as for a 2 x 2 matrix, except that the algebra on the way 
is likely to be more complicated. 


Example 6 
Find the eigenvalues and corresponding eigenvectors of the matrix 
os 1 1 
A= {4 2 1 
0 0 5 


Solution 
First, we find the eigenvalues, by finding the values of A which satisfy 
det (A — AI) = 0, that is 
2A 1 1 
1 2-A 1 = f= 0. 
0 0 5-—A 


This determinant can be evaluated in various ways. One way is to interchange 
rows 1 and 3, to obtain a form which is easier to expand: 


0) 0 5-—A 
—j} 1 2-—A l to Here we have used the 
5 \ 1 property that interchanging 
s two rows of a determinant 


changes its sign. (Unit 20, 


; ee 
Now, expanding by the top row, we obtain the characteristic equation Sction 5, Property 1.) 


—6-—A)(i- 2-4’) =0 
or (5—A)\(3-—A)\1—-—A) =0. 
So the eigenvalues are 5, 3 and 1. 


To find the eigenvectors, we solve (A — AI)x = 0, that is 


2—A 1 1 x4 
1 2-4 oS ey (5) 
0 0 5—i X3 


putting A equal to each eigenvalue in turn. 
Case 4 = 5: in this case Equation (5) becomes 
—3x, + x, +x3=0 E, 
xX, —3x,+x3=09. E, 


Note in this case that the third equation is 0 = 0. Since there are now eflectively 2 
equations in 3 unknowns there is an infinite number of solutions. 


We can find the solution to these equations, using Gaussian elimination. E, + 3E, Gaussian elimination was 
gives discussed in Unit 9. 


8 io 
—3X2 + 3x3 = 0. Ex, 


k — 
Now putting x; = k, from E,, we get x, = 5 Substituting these values into E,, we 


- get 
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k 
—3x, Me eee 


a k 
which gives x; = 7 Thus we get the solution to the equations: 


Xo =z, X54 =. 


Hence, an eigenvector corresponding to the eigenvalue 5 is 


me NR Ne 


An equally good eigenvector would be [1 1 zy. 
Case 2 = 3: putting this value of 4 into Equation (5), we obtain 
=A Xe +X, = 0 E; 
Xi —X,+ x3,=0 E, 
2x, =9 E;. 


From E3, we see that x; = 0. Substituting x; = 0 into E, and E, we see that they 
both provide us with the same information, 


xX, —x, =), 
from which we obtain a solution 
x; =k, X, = Kk, x3 = 0 
and hence [1 1 0] is an eigenvector corresponding to the eigenvalue 3. 
Case 4 = 1: finally, putting 2 = 1 into Equation (5), we obtain the equations 
X, +X, + x3 =0 


Myo xg ts eS VU 
4x, =0. 


Again, from the last equation, x3; = 0. Substituting this into the first two equations, 
we obtain the single equation 


Xi +x, =0. 
Hence 
x, =k, xX, = —k, x, =0 
is a solution to the equations, and [1 —] 0]" is an eigenvector 


corresponding to the eigenvalue 1. 


Exercise 3 
Find the eigenvalues and corresponding eigenvectors of the matrix 


1 2 6 
A=]0 + 4). 
0 0 Z 


[Solution on p. 45 | 


For a 3 x 3 matrix, the characteristic equation which we get from det (A — AI) = 0 
will be a cubic. In all the examples and exercises which we have used in this 
section, the characterstic equation factorizes easily, and has integer roots, but this 
is hardly likely to be the case in practical situations. 


In general, you will have to use some numerical method to find the eigenvalues. If 
you have access to a computer (or a programmable calculator), then the methods 

described in the rest of this uit should be used. If not, then you could always find 
the roots of the characteristic equation by using the Newton—Raphson method 


To economize on space we 
shall often print column 
vectors such as this in the 
form [4 4 1]7, using 
the transpose notation 
introduced in Unit 20. 
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described in Unit 18. Sketch the graph of the characteristic polynomial, and use 
the approximations you find from the graph as the three starting values for 
Newton—Raphson. In this way, find each of the two roots to the accuracy of your 
choice. After this, the eigenvectors can be found in much the same way as before. 


You could also use the methods described in the next section to find these 
eigenvalues—there is not a great deal of difference in the amount of work 
involved. 


In either case, this approach is only suitable if the roots of the characteristic 
equation are real. Fortunately, in practical applications, like the ones you will 
meet in Unit 24, they usually are. 


For anything larger than a 3 x 3 matrix, the whole process becomes extremely 
laborious, and you would be advised to get a computer package to help you. 


14 Nature of solutions to the eigenvalue problem 


So far in this section, all the matrices have had real eigenvalues. For instance, the 
3 x 3 matrix in Example 6 


9 1 1 
1 ,. 1 (6) 
ee 


has three real eigenvalues 5, 3 and 1. 


However, there are many matrices which do not have this property. For instance, 
the matrix in the example below has two complex eigenvalues. 


Example 7 


1 -1 
The matrix B = 4 | has eigenvalues given by 


1—A 4 


4 i 


The characteristic equation 
(1—1)*+4=0 

can be simplified (using (1 — 2)? + 4 = (1 — 1)? — (2i)”) to give 
(1 —4A+ 2i)(1 — A — 2i) = 0. 


Thus, B has complex eigenvalues 1 + 2i and 1 — 2i. 


Exercise 4 
Find the eigenvectors of the matrix in Example 7. 


[Solution on p. 45] 


Another important property which the matrices we have looked at so far have 
had, is that all the 2 x 2 matrices have had two eigenvalues, and all the 3 x 3 
matrices have had three eigenvalues. For instance, the matrix (6) above has three 
distinct eigenvalues 5, 3 and 1. More generally, it is usually the case that an n x n 
matrix has n distinct eigenvalues. However, there are many matrices which do not 
have this property, as demonstrated in the following example. 


Example 8 
5 3 
The matrix A = 3 7] has eigenvalues given by 
5—A 3 
ee 


The characteristic equation 


S -Ai-1 1) +9 20 
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can be simplified to give 
A*—444+4=0 
or (A — 2) = 0. 


In this case, we say that A has a repeated eigenvalue of 2, and that the eigenvalues 
of A are not distinct. 


When we try to calculate the eigenvectors we obtain the equation 


X14 + xX = 0. 
This equation has a solution of the form x; = k, x, = —k, and hence there is an 
eigenvector [1 1]? corresponding to the eigenvalue 2. 


Note that the matrix in Example 8 effectively has only one eigenvector 
corresponding to the repeated eigenvalue 2, whereas all the other 2 x 2 matrices 
we have met have had two distinct eigenvalues and eigenvectors. However, there 
are cases in which a repeated eigenvalue gives rise to more than one eigenvector. 
This is illustrated by the following example. 


Example 9 
By inspection, it can be seen that the eigenvalues of 
2 0 1 
A=10 2 1 
0 0 3 


are given by 
(2 —A)j?(3 — A) = 0. 
Thus, A has two eigenvalues 2 and 3. The eigenvalue 2 is a repeated eigenvalue. 


Putting A = 3 into (A — AI)x = 0, we obtain the equations 


aaa. 7 + Ag = 0 
— - x = 0. 
Hence we get an eigenvector [1 1 1]? corresponding to the eigenvalue 3. 


Putting A = 2 into (A — AI)x = 0, we obtain only the information 
Az = 0. 


So x, and x, can assume any values we choose, and these values are not related to 
each other. Thus, any vector of the form [k l 0}7 will be an eigenvector. In 
particular [1 0 Oj]? and [0 1 0)’ are both eigenvectors, and these 
are essentially distinct. More precisely we mean that the vectors are linearly 
independent. Further, the third eigenvector [1 1 1]? is essentially distinct 
from both these in the sense that all three vectors are linearly independent. This 
can be shown by forming the vectors into the columns of a matrix P, say 


1 0 1 
P=/|0 1 ioe 
0 0 1 


Then det P = 1, which is non-zero, and this is sufficient to show that the rows and 
columns of P are linearly independent. 


In general, an n x n matrix has at most n distinct eigenvalues. If there are exactly n 
there will always be n linearly independent eigenvectors, but if any of the 
eigenvalues are repeated there may or may not be n linearly independent 
eigenvectors, as illustrated by Examples 8 and 9. 


In the numerical methods we discuss later in this unit, it is assumed that the 
matrix in question has n real and distinct eigenvalues. In fact, a lot of matrices 
which occur in practice do have this property. In particular, symmetric matrices 
always have real eigenvalues and a full set of linearly independent eigenvectors. 


See Unit 20 Subsection 4.4. 


A symmetric matrix A is one 
for which A’ = A (see Unit 
20, Subsection 4.4). 
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Summary of Section 1 
The set of homogeneous linear simultaneous equations given by 
Ax = AX, (1) 


where A is a known square matrix, has non-trivial solutions for various values of 
A. These values of 4 are called the eigenvalues of A. For each eigenvalue / there 
corresponds an infinite set of solutions to Equation (1). Any non-zero solution x 
in this set is called an eigenvector. 


One way of finding the eigenvalues of A is to find the values of 4 which satisfy the 
characteristic equation 


det (A — AI) = 0. 


To find an eigenvector corresponding to an eigenvalue J, find any non-zero 
solution to the equation 


(A — AI)x = 0. 


The resulting eigenvector x is not unique. 


End of section exercise 


Exercise 5 
Find the eigenvalues and the corresponding eigenvectors for the following matrices. 


(i) & (ii) ' (iii) |! - 
Se | 4 5 i 1 1 


(iv) j|1 0 -l (v) | 3 2 2 Hint: if you get stuck, (iv) has 
1 2 1 5 9) 0 an eigenvalue 1. 
2 2 3 2 0 4 


[Solution on p. 46 | 


2 Iterative methods for finding selected 
eigenvalues (Television Section) 


2.1 Eigenvalues of matrices related to a given matrix A 


The television programme can roughly be divided into two parts. The first part of 
the programme explains the geometric significance of the eigenvalues and 
eigenvectors of a 2 x 2 matrix. The second part of the programme makes use of an 
observed geometric property of eigenvectors to construct some numerical methods 
for finding particular eigenvalues and eigenvectors. 


To understand what is going on in the second part of the programme, you will 
need to know how the eigenvalues (and eigenvectors) of A~*, and (A — pI)~* are 
related to those of A. For the purposes of this subsection, A is an n x n matrix 
with n real distinct eigenvectors, We shall call these 


ee eee ae 


Since the eigenvalues are distinct, we know from Subsection 1.4 that there must be 
a corresponding set of linearly independent eigenvectors 


%, Bas nee 

where the xs are non-zero. 

We know that the eigenvalues and eigenvectors of A are related by the equations 
AX = A553, Ax, = 12X>, ar AX, = AnXn- 

More briefly we can express all these equations as 


Ax; ce A;X; for i= fo2 4 (1) 
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We can use these equations to find the relationship between the eigenvalues of A 
and A’, as the following example demonstrates. 


Example 1 


Show that the eigenvalues of A’ are 17,13,...,42 and that the eigenvectors are the 
same as those of A. 


Solution 
We know from (1) that if A has a set of eigenvalues {1,,/2,...,4,}, then 

Ax; = A;x; fort = 2 awn (1) 
where x; is an eigenvector corresponding to 4;. 


Suppose, for each i, we left multiply both sides of (1) by A. Then we obtain 


As = AG 
= A,(Ax;) (as A; is just a number) 
= A;(A;x;) (from (1)) 
SO 
A’x; = 1?x; for i= 1,2,...,n. (2) 


This last equation tells us that A? has eigenvalues 1?,13,...,1?, with corresponding 
eigenvectors X,,X>,...,X,—-Wwhich are the same as those of A. 


So, we conclude that 
(i) A? has eigenvalues 1?,13,...,A7; 
(ii) A and A? have the same eigenvectors. 


The eigenvalues of A + gl 


Using an approach similar to that used in Example 1 we can find the eigenvalues 
of A + gI, for any number g. 


We know that 

Ax; = A;X; for i= 1,2, 8 (1) 
To obtain the eigenvalues of A + gI, we add qx; (or equivalently gIx;) to both 
sides of (1). We get 

Ax; + gIx; = 4;x; + qx; for f= 
which can be rewritten as 

(A + gI)x; = (A; + q)x; fer i= 12... , (3) 


So, we conclude that A + gI has n distinct eigenvalues 
A, + 4, Ar +4, a An + q, with corresponding eigenvectors 
X1,X2,-.-,X,——Which are the same as those of A. 


Exercise 1 


Given that A has a set of eigenvalues {/,,A3,...,A,}, state the eigenvalues of 
(i)A+51 @i)A-—2I, G@n)A-—a@L. 


[Solution on p. 47 | 


The eigenvalues of A~* 


Assuming that A is non-singular we can use the approach of Example 1 to find the 
eigenvalues of A~'. Again we start with Equation (1): 


Ax; = 1X; for i= | ee 
Multiplying by A~*, we get 
x, =/1,A~ ‘x; fort = 1,2, 2..,n 
and dividing by /,;, we obtain 
1 
A@'x,; = —x; for i =1,2,...,n. (4) 


You will need this result in 
Subsection 3.3. 


Note: A;-# 0, for otherwise x; 
(=/,A~+x;) would be zero, 
contrary to the definition of 
‘eigenvector’. 
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a 1 
From Equation (4), we conclude that A~* has n distinct eigenvalues ES aaa 

1 2 n 
with corresponding eigenvectors x,,X>,...,X,——Which are the same as those of A. 


The eigenvalues of (A — pI)~' 
In the television programme we will need to know the eigenvalues and 
eigenvectors of matrices of the form (A — pl)~' for various values of p. The 


following exercise asks you to express these eigenvalues and eigenvectors in terms 
of those of A. 


Exercise 2 


Given that A has n distinct eigenvalues 1,,/,,...,4, with corresponding eigenvectors 
X;,X2,...,X,, and that pis a number such that A — pl is a non-singular matrix, show that: 


(i) (A — pl)~? has n distinct eigenvalues 
1 : 1 
Se 


(11) the corresponding eigenvectors are the same as those of A. 


[Solution on p. 47] 


2.2 Pre-television notes 
The following facts are assumed in the television programme. 


1. The product of a matrix and a column vector gives another column vector. 
For instance, 


falll-[7] 


In the last unit we used matrix multiplication to describe a change of co-ordinate 
axes. The geometrical interpretation used in this unit is slightly different, in that 
we fix the co-ordinate axes and say that a matrix acting on a vector produces a 
new vector. The process of producing a new vector in this way is often referred to 
as a transformation. For example [3 1]? is moved to [11 7]" by the 


3 2 oe 
transformation described by the matrix 1 f as in Figure 1. In particular, the 


television programme refers to transformations. 


2. You will remember, from Section 1, that the eigenvalues of a 2 x 2 matrix A 
can be found by solving 


det (A — AI) = 0. 
As this invariably produces a quadratic equation, there are three possibilities: 
(i) A has two distinct real eigenvalues; 
(11) A has two equal real eigenvalues; 


(11) A has two complex eigenvalues. 
aS 3 Z 1 -1l 
The matrices referred to in the television programme are \ 4 and \ it 


In Section 1 you found that has 2 distinct real eigenvalues with 


eigenvectors [1 1]’ and [1 —4]’. If you have done Exercise 4 at the end of 


1 -1 
Section 1 you will have found that the second matrix 1 _ has complex 
eigenvalues. 


3. From the work in Subsection 2.1, you will need to know that if A has the set 
of eigenvalues {/,,/,,...,4,} with corresponding eigenvectors x;, x,...,X,, then 


(i) A~* has the set of eigenvalues 


1 1 1 
a 


so long as A is non-singular. 


This result is needed for the 
television programme. 


Unit 20, Section 2.4. 


i 


Figure 1 
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(ii) (A — pl)~* has the set of eigenvalues 


1 1 1 ! 
AP ae 


so long as A — pl is non-singular. 


(iii) A, A~* and (A — pI)~* all have the same eigenvectors. 


Now watch the television programme. In the second part of the programme, 
concentrate on the methods, rather than the arithmetical details. 


The following subsections all contain references to the television programme. Specific 
references will be marked with a symbol. 


2.3 Geometric interpretation of eigenvalues and eigenvectors 


(Summary of first part of the television programme) 


1 4 
vectors. In this way a vector x is found for which the vector formed by taking 


; ik 


is in the same direction as x. 


Eo 2 
1. The programme examines the effect that the matrix has on different 


In fact, it is shown that a transformation on any vector in this ‘special direction’: 
(i) leaves the direction unchanged 


(ii) stretches the vector to five times its original length, that is 


3 ys 
k = IM (See Figure 2.) 
1 4 
2. Other matrices are examined to see if they have associated ‘special directions’. 
oe. Sl 
It is found that some do, and some do not. For intance, the matrix 1 1 has 
no associated ‘special direction’. It is, in fact, a multiple of the rotation matrix 
co a sin 
ee ‘iii 5 es eens 
Je ee 4 4 
1 le n | 
sin — COs — 


A Jue 


: 7 ; 
and has the effect of rotating every vector it operates on through an angle 1 anti- 
clockwise, and enlarging it. 


3. To discover whether it is possible to predict which matrices have ‘special 
directions’ and, if they do, by how much the vectors in this direction are stretched, 
the programme looks at the algebra of what has been happening. We have been 
looking for ‘special vectors’ x such that 


Ax = 1x, 


where 4 is the stretch factor. This equation is just the algebraic eigenvalue problem 
discussed in the first section. 


1 1 
are complex. For a geometric interpretation, however, we need a real eigenvector. 
As this matrix has no real eigenvectors, it has no special direction either. 


; : 1 -l 
Exercise 5 of Section 1 shows that the eigenvalues and eigenvectors of | 


— 3 : 
In the pre-television notes we recalled that 1 i has two eigenvectors: 


[1 1]? and [1 —+4]". The television programme checks geometrically that 
[1 —4]? is a special direction, and that the enlargement in this direction is 
indeed the corresponding eigenvalue. 


TV 21 


Figure 2 
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The conclusion we come to is that for some 2 x 2 matrices A, there are ‘special 
vectors’ x, which when transformed do not change their direction. These special 
vectors, when they exist, are eigenvectors of A. The factor by which each of these is 
stretched can be found by calculating the equivalent eigenvalue. 


So, if A has real eigenvalues 4 with corresponding eigenvectors x, then 
(i) x will be in a ‘special direction’ 


(ii) 2 will be the amount by which x is stretched. 


24 Direct iteration 


oO 
The television programme takes a matrix A which is known to have two real ; 
eigenvectors. Arbitrary vectors y are taken, and transformed by A to give new 
vectors Ay. All appear to be closer than the corresponding vector y to one of the 
‘special directions’ (see Figure 3). This holds provided that y was not in a special 
direction to start with, for then Ay is known to be in the same direction as y. To 
see why the new vector is closer to the special direction than the old one let us 


3 2 
look at the matrix A = 1 4/ which we know has eigenvalues 5 and 2 with 
corresponding eigenvectors [1 1}’ and [1 —4]’. 


special direction 1: 
stretch factor 5 


special directions 


Figure 3 
= 
fen 
special direction 2: 
stretch factor 2 
Zz, = 
Figure 4 
Any two-dimensional vector y not in one of the special directions can be expressed 
as the sum of vectors in the two special directions. As vectors in these directions 
are eigenvectors, y can simply be expressed as the sum of two eigenvectors of 
suitable lengths: 
y =X, + X>. (5) | Remember, eigenvectors are | 
Seas. not unique. 
Multiplying (5) by A, we find that 
Ay = Ax, + Ax, Ax, = 5x; and Ax, = 2x> as 


x, and x, are eigenvectors. 
= 5x, + 2x. 1 2 


Now the stretching in the x, direction is greater than that in the x, direction. So 
Ay is closer than the original vector y to the x, direction. 


More generally we can consider any 2 x 2 matrix with distinct real eigenvalues /, 
and /,. For then an arbitrary vector y (not itself an eigenvector) can be expressed 
as the sum of two eigenvectors: 


y =X, + Xp. (5) 
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Hence 
Ay = Ax, + Ax, 


= A,X, + Ax, (since x, and x, are eigenvectors). 


Now the vector Ay is closer than y to one of the special directions specified by x, 
and x,. Which one it is closer to depends on which of 1, and 1, has the larger 
magnitude. In general, Ay will be closer than y to the eigenvector which 
corresponds to the eigenvalue of largest modulus. 


Provided that the eigenvalues are real, there are only two exceptions to this 
statement: 


This vector is referred to in 
the television programme as 
the dominant eigenvector. The 
corresponding eigenvalue is 
called the dominant eigenvalue. 


(i) if y happens to be an eigenvector, then Ay must be in the same direction as y 


(11) if A has two eigenvalues of equal modulus, then the statement fails because 
there is no dominant eigenvector. 


However, the fact that, in general, Ay is closer than y to the dominant eigenvector 
means that if we multiply the vector Ay by A again, the new vector formed will be 
even nearer to the dominant eigenvector. This observation gives rise to an iterative 
method for finding the dominant eigenvector. If we start with a vector y, (not an 
eigenvector), and form the sequence of vectors 


Yi =Ayo 
Yo =Ay; 
Yr+ | eee Ay, 


then each successive vector will be nearer than the previous one to an eigenvector 
corresponding to the eigenvalue of largest modulus. If we continue in this way 
then, for large enough 7, y, will effectively be an eigenvector of A. 


This process is called direct iteration, and is demonstrated in the following 
example. 


Example 2 


3 4/ and an arbitrary non-zero vector yo = [1 ie 


(not an eigenvector of the matrix), then the iterative scheme above can be used to 
find an estimate for the dominant eigenvector. 


BOURKE 
IEE 
EIB 


and so on. 


3 
Given a matrix A = 


To see how quickly the process is converging, we can check the ratios of the 
elements. We can do this by dividing by the element of largest modulus. Thus 
dividing y, by 43 gives y, = 43 [0.67 1]’. Similarly, dividing y; by 259 gives 
y3 = 259 [0.67 1]", so, it seems likely that a required eigenvector of A is 
approximately [0.67 ry. 


It is also possible to find the dominant eigenvalue using this process—but first we 
shall amend the process into a more usable form. You may have noticed that the 
numbers in the example above were getting larger and larger. This is more than a 
mere inconvenience. If the numbers were allowed to go on growing in this way, it 
is possible for the numbers in the successive column vectors y, to get too large 
even for a computer. However, as it is only the ratio of the numbers in the vector 
that matter, we can scale each vector at each step so that the elements in the 
vectors stay a reasonable size. There are many ways of scaling, but the one usually 
used is to divide by the element of largest modulus at each stage of the iteration. 
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For instance, if at some stage the iteration gave the vector y’'= [2 —1]’ we 
would divide by 2 to obtain y = [1 —4]’. Similarly, if the iteration gave the 
vector z’= [1 —3]? we would divide by —3 to obtain z= [—4 1]. In this 
way one of the elements in the scaled vector will be 1 and none of the other 
elements will be greater than 1. This is demonstrated geometrically in the 
television programme, by showing that each 2-dimensional vector formed is in turn 
scaled onto a point on the square (see Figure 5). 


Figure 5 


The following example demonstrates this scaling process. The notation y, is 
adopted for the unscaled vector, and y, for the vector after it has been scaled. 


Example 3 


Re-working the last example, but this time dividing through by the element of 
largest modulus at each stage, we get a new sequence of vectors: 


3 2111 5 
= | = | All results were worked to full 
3 4} 1 7 : accuracy, but are recorded 


Dividing by 7 gives y, = (0.714 | i}. here to three decimal places 


for convenience. 
ee ar 3 2 || 0.714 S 4.143 
nl TE ee 6143 | 


Dividing by 6.143 gives y, = [0.674 1). 


iy ya [3 21/0974] _ [4023 
ee 
Dividing by 6.023 gives y; = [0.668 Ey . 


This time, the numbers have stayed under control. To see that [0.67 i} ae 2 
reasonable approximation for the dominant eigenvector, look at y, and y3; these 
are both [0.67 1] to two significant figures. 


In general, to see if a required degree of convergence has been reached, we check 
that each element of y,,, — y, is small enough in magnitude to satisfy the 
accuracy we require. For instance, in the example above 

¥3 — y2 = [—0.006 0)", so the largest magnitude is 0.006 which is slightly 
larger than the value 0.005 required for two significant figure accuracy. In practice 
therefore, we really should do another iteration—and check that the largest 
magnitude in y, — y; is less than 0.005 before we can be satisfied that the given 
answer is indeed correct to two significant figures. The dominant eigenvector is 
actually [4 1. 
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This process also gives us an approximation for the dominant eigenvalue. From 
Part (iii) of the iteration, we know that 
y3 _ AY2 


13 6023 60S 


and, since y, ~ y3, we have 


AYy3 
6.023 


y3=> 


This can be re-written as 
Ay; ~ 6.023y3, 


which tells us that as well as y; being an approximation for the dominant 
eigenvector, 6.023 must be an approximation for the dominant eigenvalue. 


In this case, three iterations have produced a good approximation. 


Thus we obtain the following iterative scheme. 


Procedure 2.4: Direct iteration 


To find the eigenvalue of largest modulus and corresponding 
eigenvector of a given matrix A, start with an arbitrary non-zero 
vector Yo (not an eigenvector) and then 


1. form Yea a Ay,, 


/ 

Yr+1 
2. form y,4, =——, 

r+1 


where «,,, is the element of largest modulus in y, +. 


B64 and 2 for r = 0, 1, 2.4% 


Then, for sufficiently large r, «,,, will be a good approximation 
to the eigenvalue of largest modulus of A, and a corresponding 
eigenvector will be approximately equal to y,+1. 


This scheme will work, so long as 
(i) the eigenvalues are real, 
(ii) the moduli of the eigenvalues are distinct, and 


(iii) the initial vector yo is not an eigenvector. 


Exercise 3 
Use the procedure above to find y,, y2 and y; for the matrix 


— 4 Z 
So Foe 
starting with yo = [1 Er 


Hence, estimate the eigenvalue of largest modulus, and a corresponding eigenvector. 


[Solution on p.47 | 


This procedure for finding the eigenvalue of largest modulus is only demonstrated 
here for a 2 x 2 matrix. However, it works equally well on larger matrices which 
have real distinct eigenvalues. The advantage of this method over the approach in 
Section 1 is that iterative schemes like this one are suitable for use on a computer, 
whereas methods based on the characteristic equation are not. 


You will have the opportunity to use the method on larger matrices when you use 
the computer package described in Section 4. In the meantime, you can see the 
computer printout we got using the computer package for the 3 x 3 matrix 
discussed in the television programme, at the end of this section. 


This is often referred to as the 
power method. 
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To discover something about the rate of convergence of direct iteration, we look at 
the algebra of the crude method described initially where scaling was not used. 


Let A be an arbitrary matrix with eigenvalues satisfying the conditions following 
Procedure 2.4. We start by writing the initial vector yy in terms of eigenvectors of 
A: 


Yo = X1 + Xp. 
The first iteration gives 
Yi = Ayo 
om A(x, 4 X>) 
— AiX4 ss A>X>. 


Similarly, the second iteration gives 
Y2 = A(A,xX, + A2X2) 
= /,(Ax,) + 42(Ax2) 
= J2x, + 13x). 
In general after r iterations we obtain 
y, = AiX, + AbX>. 


If we rewrite this in the form 


Ao\" 
in 


and assume that /, is the eigenvalue of largest modulus, then |/2/A,| < 1, and so 
(A,/A,)’ will become very small as r gets larger. Thus for large enough r, 


i= M x, * 


yr = AX}. 


So y, is approximately an eigenvector corresponding to the eigenvalue of largest 
modulus. 


The working above shows that if the eigenvalues are well separated, then |A2/A,| 
will be much less than 1, and hence convergence will be rapid. 


For instance, in Example 3, you saw that we obtained a good approximation for 
both the eigenvalue of greatest modulus and its corresponding eigenvector in three 


3 2 
iterations. This was because the eigenvalues of 3 7 are 1 and 6 so (1/6)’ gets 


small very rapidly. Unfortunately this is not a fact we know in advance, as the 
eigenvalues are unknown. 


The algebra above gives some indication as to how we would examine 
convergence, and derive the method for larger matrices, where the geometry 
discussed in the television programme is inapplicable. 


2.5 Inverse iteration 


The direct iteration method (or power method) described in the last subsection has 
the obvious disadvantage that it only finds one eigenvalue (and eigenvector)—the 
one of largest modulus. 


However, in this subsection, we modify this method, and use the modified method 
to find other selected eigenvalues of a-matrix A. First, let us consider the iterative 
scheme : 


Y+1=A'y, (6) 


where A is some non-singular matrix, and yy is some arbitrary non-zero vector 
which is not an eigenvector of A. 


From the work in the previous subsection, we know that this scheme will produce 
a sequence of vectors y,,y2,... which will tend to the dominant eigenvector of 
A~'. Further, if we scale this scheme, as in the previous subsection, we will obtain 
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the eigenvalue of A~* with largest modulus. We know from Subsection 2.1 that 
the eigenvalues of A and A~* are related, so the fact that we can find an 
eigenvalue p of A~’, means that we know 1/y will be an eigenvalue of A. Now, we 
know that p is the eigenvalue of A~' with largest modulus, so 1/u must be the 
eigenvalue of A with smallest modulus. 


This modified scheme is called inverse iteration and is implemented, including 
scaling, in the following example. 


Example 4 
We shall use the scheme described above to find the smallest eigenvalue of 


Solution 
The inverse of A is 


if 4 -2 
ee Sees 
: = _; ) 


Choosing an initial vector yp = [1 1]7 and using scaling at each stage, we get 


- If 4 re 
0) or i l-| eet 


Dividing by 0.333 we obtain y, = [1 OT. 


i ea 
a) ind 3 3 e gs 


Dividing by 0.667 we obtain y, = [1 —0.750]". 


eee =o 1-4 
wm w3 63 3 O75 ee 
Dividing by 0.917 we obtain y; = [1 —0.955]’. 


Thus, an approximation to the eigenvalue of A~* with largest modulus is 0.917 
and a corresponding eigenvector is [1 —0.955]’. Hence an approximation for 
the eigenvalue of A with smallest modulus is 1/0.917 = 1.091, with a corresponding 
eigenvector [1 —0.955]", for we know from Subsection 2.1 that A and A * have 
the same eigenvectors. The actual eigenvalue is 1, with a corresponding eigenvector 
[1 —1]’,so this approximation after three iterations is not bad. 


We summarize this first version of the inverse iteration procedure below. 


Inverse iteration (first version) 


Given a non-singular matrix A, and an arbitrary non-zero vector 
Yo (not one of the eigenvectors of A) 


t. find ¥!.; =A ‘yz See practical note after next 


: exercise. 
Yr+1 

2. Wad y,.¢ = : 
Op +1 

where «,,, is the element of largest modulus in y,,. 


Do 1 and 2 forr = 01. 2..% 


1 

Then for sufficiently large r, will be a good approximation 
r+t 

for the eigenvalue of A with least modulus, and the corresponding 


eigenvector will be approximately equal to y,4+ 1. 
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Exercise 4 


7 3 
Use the procedure above with A = : a and yy = [1 1]" to obtain y,, y, and y3. 


Hence find an approximation to the eigenvalue of A with least modulus and a 
corresponding eigenvector. 


[Solution on p. 47 | 


Practical Note 


For matrices larger than 2 x 2, the calculation of the inverse matrix is fairly 
laborious, and the step y,,; = A “y, is done by solving the equation Ay,+, = y,. 
However, as the inverse of a 2 x 2 matrix is very easy to calculate, we advise you 
to use A~* when doing these examples by hand. 


A more versatile form of inverse iteration 


The scheme described above can be generalized to obtain a procedure which will 
find the eigenvalue of A which is closest to some chosen value p. 


The crude scheme which we use is 
Yr+1 = (A — pl)“’y,. (7) 


This will, with scaling, certainly determine the eigenvalue y of (A — pl)~! which 
has largest modulus, and this enables us to determine an eigenvalue of A. We 


1 
know from Subsection 2.1 that if 2 is an eigenvalue of A, then is an 
Bae 


eigenvalue of (A — pI)~'. So we can obtain this eigenvalue 1 from 


Now, as pu was the eigenvalue of (A — pI)~* with largest modulus, 1/u must be the 
eigenvalue of A — pI with smallest modulus. In other words 1/ is the eigenvalue 


1 
of A — pI closest to zero and so — + p is the eigenvalue of A closest to p. The 
Ll 


following example demonstrates the use of the Scheme (7), with scaling at each 
stage. 


Example 5 


3 Zhe 
Suppose A = 3 a Find the eigenvalue (and corresponding eigenvector) 


nearest to 7. 


Solution 


Putting p = 7 gives 


So we shall implement the Scheme (7), scaling at each stage as before, with an 
initial vector yp = [1 i}. 


. 2 ft -093 
w) vie -7 ef eee 


Dividing by —1.167 we obtain y, = [0.714 1]°. 
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: A a ee 
m) = -73 mt ia “| tal 


Dividing by — 1.024 we obtain y, = [0.674 ig a 


=. Ip en ee 
Oe a 4 a 


Dividing by — 1.004 to obtain y3 = [0.668 aS 


Thus, an estimate for the eigenvalue of (A — pI)~* with largest modulus is — 1.004, 
and [0.668 1}’ is the corresponding eigenvector. Hence, an estimate for the 


1 
A ‘2. 
eigenvalue of A nearest to 7 is “T1004 


is [0.668 1]’ for we know that A and (A — pI)! have the same eigenvectors. 


Again, these turn out to be good estimates for the actual eigenvalue 6 and 
eigenvector [2/3 . 

This scheme, as formally set out below, is the one we shall refer to as inverse 
iteration. 


Procedure 2.5: Inverse iteration 

To find the eigenvalue closest to p and the corresponding 
eigenvector of a given matrix A, start with an arbitrary non-zero 
vector yp (not one of the eigenvectors of A) and then 


1. find Ve+1 = (A ar pl)~'y,, 
Yr+1 


Z. find vrei os 5 


br+4 


where «4, is the element of largest modulus in y,, ;. 
01 ane 2 tor f= 01.2... 


Then for sufficiently large r, + p will be a good 
ae 


approximation for the eigenvalue of A nearest to p, and y,., will 
be a good approximation for the corresponding eigenvector. 


This result in fact covers the case of the previous boxed result as, by putting p = 0, 


we get the eigenvalue nearest to 0, that is, the eigenvalue of A with smallest 
modulus. 


Exercise 5 


; 2 1 
Find y;, y2 and y3 using the scheme above, where A = ; , p = 2, and using a 


starting vector yo = [1 0]’. Hence, estimate the eigenvalue of A nearest to its 
corresponding eigenvector. 


[Solution on p. 48 | 


A note about choosing p for a 3 x 3 matrix 


For a 3 x 3 matrix the largest and smallest eigenvalues can be found by iterating 
with A and A~* respectively. The remaining eigenvalue can be found by iterating 


with (A — pI)~* provided a suitable value can be chosen for p. 


+ 7 = 6.004. The corresponding eigenvector 


See practical note on page 23. 
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In the television programme the three eigenvalues of the matrix 


2d 
0.5 
Be 


25 


are found, using the computer package EIGSOL described in Section 4. The 
eigenvalues of largest and smallest modulus were found to be approximately 9.91 
and — 6.24 respectively. 


-9.9] 


Figure 6 


eigenvalue 


-6.24 


6.24 


eigenvalue 


\ 
9.91 


Having found these two eigenvalues the choice of a suitable p can be cut down to 
two regions (see Figure 6). The remaining eigenvalue cannot lie between — 6.24 
and 6.24—for otherwise it (and not 6.24) would be the eigenvalue of smallest 
modulus. Nor can it be larger than 9.91 or smaller than —9.91, for then it (and 
not 9.91) would be the eigenvalue of largest modulus. So the remaining eigenvalue 
must either lie between 6.24 and 9.91, or between — 6.24 and —9.91. So try 

p = (6.24 + 9.91)/2 ~ 8, or if this does not work try p = —8. One of these values 
of p will give you the remaining eigenvalue. 


In general, if 2; and A, are the eigenvalues of largest and smallest modulus 


respectively for some 3 x 3 matrix, then choosing p to be either 


— [Waal + al 
Ps 


~—TAs| 


Figure 7 


—{A2| 


0 


dil + Vel 
2 


| is bound to produce the remaining eigenvalue (see Figure 7). 


Finally, here is the computer output we got when finding the three eigenvalues 


(and eigenvectors) of the matrix used in the television programme. 


(i) 


1 

1 
0.98047 
0.95930 
0.95542 
0.95439 
0.95415 
0.95410 
0.95408 
0.95408 


DIRECT ITERATION 


Yn 


—0.00848 


0.38486 
0.32369 
0.32176 
0.31986 
G.31955 
0.31947 
0.31944 
0.31944 
0.31944 


Eigenvalue of largest modulus is 9.90552 
Found in 45 iterations 


Corresponding eigenvector 


0.95408 
0.31944 
1 


0.52542 
0.97957 
1 


ew eg est =f aft ott = 


ALPHA 


11.8000 
ee 
9.87006 
9.91254 
9.90566 
9.90570 
3.90555 
9.90553 
9.90552 
9.90552 
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(ii) INVERSE ITERATION (WITH p =0) 

n Yn 

1 1 — 0.31642 0.22687 
11 -0.41539 1 0.42106 
21 —0.31521 1 0.40413 
31 —0.29990 1 0.40010 
41 -0.29728 1 0.39940 
51 —0.29683 1 0.39928 
61 -0.29675 1 0.39926 
rg -—0.29673 1 0.39925 
81 -0.29673 1 0.39925 
86 -0.29673 1 0.39925 


The eigenvalue of maximum modulus of the inverse of 


(A —0.00 #1) is —0.16037 
Hence the eigenvalue nearest to 0.00 is —6.23567 


Found in 86 iterations 


Corresponding eigenvector 


—0.296735 
1 
0.39925 
(iii) INVERSE ITERATION (WITH p =8) 
n Yn 
1 —0.50743 -—0.28218 1 
2 1 0.35244 0.02694 
3 1 0.35402 -0.22438 
4 1 0.35261 —U1SToU 
5 1 0.35298 =O:16222 
6 1 0.35286 —0.15472 
¢ 1 0.35290 —0.15696 
8 1 0.35289 -—0.15629 
9 1 0.35289 —0.15649 
10 1 0.35289 -—0.15643 
11 1 0.35289 -—0.15645 
12 1 0.35289 —0.15644 
13 1 0.35289 —0.15644 
14 1 0.35289 —0.15644 


The eigenvalue of maximum modulus of the inverse of 
(A — 8.00 #l) is —1.75484 
Hence the eigenvalue nearest to 8.00 is 7.43015 


Found in 14 iterations 
Corresponding eigenvector 


1 
0.35289 
—0.15644 


Summary of Section 2 
1. If A is a matrix with eigenvalues 4,,/,,...,4, then: 
(i) A+ ql has eigenvalues 4; + g, 42 + q,...,4n + G3 


1 1 
(ii) A~* has eigenvalues [re Ss so long as A is non-singular ; 
1 2 n 
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ALPHA 


0.16059 
—0.14796 
—0.15832 
—0.16001 
=0.1G031 
—0.16036 
—0.16037 
—0.16037 
—0.16037 
—0.16037 


ALPHA 


0.52271 

2.72664 
—1.41021 
—1.88353 
= 1.71900 
—1.76578 
—4,76159 
—1. 15082 
—1.75455 
—1.75493 
— 1.75482 
—1.75485 
—1.75484 
—1.75484 
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1 1 
Ay —p A, ae ee ad 


(iii) (A — pl)~* has eigenvalues so long as A — pl is 


non-singular; 


(iv) A, A + gl, A~‘ and (A — pI)~! all have the same eigenvectors. 


2. Ifa2 x 2 matrix A has two real eigenvectors, then there are special directions 

x such that the vector Ax is in the same direction as x. These special directions are 
given by the eigenvectors of A. The factor by which A stretches the original vector 

X is given by the corresponding eigenvalues of A. 


3. The procedures summarized below determine single eigenvalues of a matrix A. 
In all cases we start with an arbitrary non-zero vector y, (not an eigenvector). 
“+1 18 the element of largest modulus in U4 


Iterative scheme to be used 
for 7 = 04.9 


Procedure Result for large r 


Direct Form y,,, = Ay, %,41 approximates the 
iteration y’ eigenvalue of largest modulus. 
rt 
aa set Yi = Geax Y-+1 approximates the 
corresponding eigenvector. 
Inverse Form y;+, =A~'y, approximates the 
iteration ; Or 1 
: ¥r+i eigenvalue of smallest modulus. 
and set y,,, = 
ee Y-+1 approximates the 
corresponding eigenvector. 
Modified Form y,,, = (A — pl)“ 'y, + p approximates the nearest 
inverse wre ae h I 
‘ieranian ees eigenvalue to a chosen value p. 


Y-+1 approximates the 
corresponding eigenvector. 


Note that the second iterative scheme is a special case of the third scheme with 
p=, 


The above schemes will only work if the eigenvalue we are looking for is real and 
distinct from all the others. 


4. Fora 3 x 3 matrix, if the eigenvalues of largest and smallest modulus are Ay 


— Ay| + |A 
and 4,, then to find the remaining eigenvalue, choose p = + | ll 
End of section exercise 
Exercise 6 
For the matrix 
— 39 40 
- | — 20 = 
use direct and inverse iteration with initial vector Yo = [1 0]}' to find: 


(i) the eigenvalue of A with largest modulus; 
(ii) the eigenvalue of A nearest to 2. 


In each case, stop the method when each element in two successive approximations y,,, 
and y, agree to two significant figures. 


[Solution on p. 48] 
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3 Decomposition methods for finding all 
the eigenvalues 


3.0 Introduction 


In Section 2, you saw the development of a class of iterative methods, which can 
be used to find selected eigenvalues (and associated eigenvectors) of a given matrix 
A. Often, however, in practical problems, all the eigenvalues are required. 
Although it would be possible in theory to find these eigenvalues by repeated use 
of inverse iteration with suitable values of p, it is in fact more practical to use a 
method which is specifically designed to find all the eigenvalues simultaneously. 
Most methods of this type currently in use depend on suitable ‘decompositions’ in 
which the given matrix A is expressed in the form PQ for suitable matrices P and 
Q. It is for this reason that Subsection 3.1 concentrates on matrix decomposition. 
Subsection 3.2 deals with some algebraic results which will help you understand 
the method described in Subsection 3.3. This is a method which can be used to 
find all the eigenvalues of a matrix simultaneously. 


The exercises that you will be asked to do in this section are mostly on 2 x 2 
matrices. This is simply to cut down the amount of hard computation you are 
required to do. Rest assured that the processes described will work equally well for 
larger matrices, and you will have the opportunity to do some of these using the 
computer package in Section 4. 


3.1 Matrix decomposition 


Decomposition of matrices is analogous to factorization of numbers. For instance 
much in the same way that 12 has a factorization 12 = 3 x 4 you can check that 


1 1 
the matrix \ | has a decomposition 


6 


} lst, 0 I, 1 (1) 
1 6 1 ee 0 es 

Now, some numbers have more than one factorization. For instance 12 = 3 x 4, 
or 12 = 6 x 2. Similarly it is possible for a matrix to have more than one 


decomposition. For instance, you can check that the matrix above has another 
decomposition 


it l-lt eos @ 


In this subsection, we discuss just one possible decomposition that a given matrix 
may have. In certain circumstances, a square matrix A can be decomposed as the 
product 


A = LU, 


where L is a lower triangular matrix, which is chosen to have ones along its main 
diagonal, and U is an ordinary upper triangular matrix. This is often referred to as 
an LU decomposition, and an example is given below. 


Example 1 
4 1 1 a oe 
<5 5 ii= 6 
ne eee 0™ 
. = U 


Another example of an LU decomposition is the 2 x 2 decomposition (2) above. 


A simple approach to finding an LU decomposition is just to assume it exists! The 
following example demonstrates this technique. 
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Example 2 
To find the LU decomposition of 


5 3 
Lc 
bs 


assume that 


A = L U 


To find the unknown elements a, b, c and d we shall just multiply L and U 
together, to obtain 


2 k b c 
; 4 = e ac + iI ©) 


Then equating these two matrices, element by element, we get from the first row 
b= 5. C3, 

and from the second row 
abp-= 3: ac+d=5. 


Hence, a = 2 and d = +8, and the required decomposition is 


Eb sl-L illo a 


Exercise 1 
Find the LU decomposition of the matrix 


6 2 
A= : 
b 2 
[Solution on p. 48 | 
Not all matrices have an LU decomposition. Given a particular matrix, which has 


no LU decomposition, the arithmetic will break down when we attempt to do the 
decomposition, as demonstrated in the following example. 


Example 3 


1 
Show that ‘| does not have an LU decomposition. 


Solution 
If there were an LU decomposition we would have 


0 1 b ¢c 
= f 3 
e ‘| i ac + 7 fe Bhove) 
Equating coefficients, from the top row we obtain b = 0 and c = 1. Equating 
coefficients in the second row, we obtain ab = —3 and ac + d = 5. Now, it is not 


possible for ab to be —3, as b is zero. Thus, the arithmetic process breaks down, 
and we know that no such decomposition exists. 


In general, a 2 x 2 matrix has an LU decomposition 


A114 
Remember: 


det A = @11822 — @1724>}. 
a21 


so long as aj, is not zero. 


You can check the result above by multiplying out the right-hand side. This result 
has been boxed and labelled (4) so that we can use it in later examples. 
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To find an LU decomposition of an n x n matrix 


Qy4 oe Ain 1 Uy — Urn 


=| ~ a Notation: all the matrix 


0 elements replacéd by Q are 
Gn ae Qnn lt ln2 Ce 1 Unn Zero. 


is in principle no more difficult than for a 2 x 2 matrix. It just involves more 
arithmetic. There are many computer packages that will do this for you. 


Exercise 2 
Find the LU decomposition of the matrix 
1 2 3 
A=tt 6 4). 
2 1 4 


[Solution on p. 48 | 


3.2 Some algebraic results 


To gain some insight as to why the process described in the next subsection works, 
we shall need some properties of eigenvalues. For easy reference, these will be 
referred to as theorems. 


The first theorem concerns the eigenvalues of the matrices AB and BA. The 
theorem asserts that these two matrices have the same eigenvalues, as the following 
example illustrates. 


Example 4 


1 0 2 1 
Let A = I < and B = k f Then 
1 01] 2 1 2 1 
= k F 4 = ‘ 5 


So the eigenvalues of AB can be found by solving 


2—A 1 


2 5-7 


That is, 
(2—A)3—A)—2=0 (5) 


or 47 —5144+4=0 
Le. (A —4)A — 1) = 0. 
Thus the eigenvalues of AB are 4 and 1. 
Now consider BA. 
2 ii 0 3 1 
BA = ae 
iS ht 

and so the eigenvalues of BA can be found by solving 
3-—A 1 

fe 2-—A 


This gives (3 — 4)(2 — 4) — 2 = 0, which is the same as Equation (5). So the 
eigenvalues of AB are also 4 and 1. 


=o. 


Hence, in this case, AB and BA have the same eigenvalues. In general, for any 
matrices A and B, the eigenvalues of AB and BA are the same, but their 
eigenvectors are usually different. 
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Exercise 3 
Show that AB and BA have the same eigenvalues if 


1 0 1 4 
a=; and B=| j “at 


[Solution on p. 49 | 


The general result can be obtained from the definition of an eigenvalue. Let A and 
B be any square matrices. If J is an eigenvalue of AB with corresponding 
eigenvector x then 


ABx = 2x. 

If we multiply this equation on the left by B, we get 
BABx = /Bx. 

Now, Bx is just a column vector. So letting y = Bx, we obtain the equation 
BAy = jy. 


From this equation, we can see that 2 is also an eigenvalue of BA with a 
corresponding eigenvector y = Bx. A similar argument, starting with BA, shows 
that if A is an eigenvalue of BA then 2 is also an eigenvalue of AB. Thus we have 
derived the following theorem. 


Theorem 1 


For any square matrices A and B. 


(i) AB and BA have the same eigenvalues: 


(il) If x is an eigenvector of AB then Bx is an eigenvector of BA. 


This result sometimes appears in other forms. For instance, suppose P and X are 
Square matrices and that P is non-singular. If we replace A by XP and B by P“! 
in the result above, then AB = (XP)P~! = X and BA = P~'XP. Thus we can 
derive the following form of the result above. 


Theorem 2 This result will be used in 


For any square matrices P and X with P non-singular: Unit 22. 


(i) X and P~'XP have the same eigenvalues: 


(ii) If x is an eigenvector of X then P~'x is an eigenvector of 
| aD, ¢ 


Some more results which we shall need are much easier to show. 


Theorem 3 
The eigenvalues of the diagonal matrix 


wie ) 


A= os 


v 


ATC G1 1,4 75,...,Anp. 


This follows immediately from finding the eigenvalues from det (A — AI) = 0. 
Q11— A @) 


Ann — A 


— (a1, —-A)(az2 ‘id cae (Ann ey 4). 


det (A — AI) = 
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Hence we find that the eigenvalues of A are 4,1,422,.-.,@nn- In a very similar 
fashion we can derive the following more general result. 


Theorem 4 
The eigenvalues of the upper triangular matrix 


Qin 


Ann 


are GQ; 14,422,---54nn- 


This again follows from putting det (A — AI) = 0. Now the matrix 


ay,—4 ao Gin 
A-AL= 0 : 
Ann — A 
is upper triangular and so its determinant is equal to the product of the diagonal Determinants of upper 
elements. That is triangular matrices were 
discussed in Unit 20, 

det (A — AL) = (a4, — A)(Q22 — A). - (Gun — 4). Subsection 5.4. 
Hence, once again we obtain the result that the eigenvalues of an upper triangular 
matrix are its diagonal elements a; 1,422,...,Qnn- 
Exercise 4 


State the eigenvalues of the following matrices 


(i) be . Ce © ee ee | 
ee 0 sa 2s 


0 0 9 
[Solution on p. 49] 


The last theorem we state in this subsection is a result you will require for Unit 22. 


Theorem 5 

Suppose that a square matrix A has n distinct eigenvalues 
A1,A3,+++,A,. Let P be a matrix whose columns are corresponding 
eigenvectors X;,X2,...,X, Of A. Then 


Ay 
P-'AP = 4a 0 


dr 


If the eigenvalues A,,A2,..,/, are not all distinct, the result still 
holds provided we can choose a corresponding set of 
eigenvectors which are linearly independent. 


This result is illustrated in the following example. 


Example 5 


5 Z 
Using the results of Example 5 in Subsection 1.2 we know that A = 5 - has 
eigenvalues 3 and 7, with corresponding eigenvectors [—1 1]? and [1 Pfs 


So, we know that 


2 sabe slide 
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We can put these two equations as a single matrix equation in the following way: 


peeeeee = jj geeeeee 


Ee eesti 
| 


So letting P be the matrix whose columns are the eigenvectors of A, we have that 


Sai 0 
AP =P 
7} 


3 0 
PAP = ; 
. a J 


This matrix just has the eigenvalues of A on the main diagonal. 


Or 


The general result of Theorem 5 will not be proved in this unit, but the proof is 
very similar to the example above. 


Exercise 5 


2 | 
Given that A = : find matrices P and D such that 


PAP =D, 
where D is a diagonal matrix consisting of the eigenvalues of A. 


[Solution on p. 49 ] 


3.3. The LR method 


The LU decomposition of a matrix is used in a process called the LR algorithm 
for finding eigenvalues. First, suppose we have a matrix Ag, and we find its LU 
decomposition. 


Ay = LoUo. 
If we then form a new matrix 
A; = UoLo, 
we know that 
(i) In general Ay and A, will be different matrices, as LpQUy 4 UoLo. 
(ii) From Theorem 1, Ay and A, have the same eigenvalues. 


The next example illustrates how repeated use of this process leads to a method 
for approximating the eigenvalues of a matrix. 


Example 6 
Suppose we want to find the eigenvalues of the matrix 


3 
A= : 


We start by finding the decomposition Ay = LoU : 


1 01} 5 3 
Ay = 06 4 0 3 “i (from Equation (4) of Subsection 3.2) 


and then form a new matrix A, = UoLo: 


ee eee 
ai=|) mM ‘ 


_ | 68 3 
- 19a 457 
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We know from Theorem 1 that Ay and A, have the same eigenvalues. By 
repeating this process we can form a sequence of matrices A,,A3,..., each of 
which has the same eigenvalues as Ay. Let us see what happens when we do this. 


The LU decomposition of A, is 


1 0 ][ 6.8 3 
i= ies 
= ae all 0 se 


(from Equation (4) of Subsection 3.2), 


and so the next matrix A, is given by 


6.8 3 1 0 
a 0 an Be N 


_ [7.647 3 
(0.664 2.353 | 


The next iteration produces the LU decomposition 


1 0 ][ 7.647 3 
ee ee All 0 aa 


and so Az; is given by 


7.647 3 1 0 
ea 0 a ae ‘ 

_ | 7.908 3 

~FOi1s2 T0027 
It is worth noting that in this sequence of matrices A,,A>,A3,..., which we are 
forming, the bottom left-hand element (set in bold type) appears to be getting 
smaller. If this continues to happen in each subsequent iteration, then sooner or 
later, some matrix A, which we form will virtually be in upper triangular form. If 
this happens, then we know from Theorem 4 that the eigenvalues of A,, and hence 


of Av, will be the values on the diagonal. Performing one more iteration, we 
obtain the decomposition 


oe 1 0 117.908 3 
So 1 0 2.023 |’ 
ae 


and so A, is given by 


7.908 3 1 0 
ee 0 — pr 1 
We Be oa: S 
4 0Ga7 . 2023 | 
Again you can see the bottom left-hand element getting smaller. If this element is 
taken to be approximately zero, then the matrix is approximately diagonal, and so 


by Theorem 4 the eigenvalues are approximately equal to the elements on the 
main diagonal. 


Thus we obtain approximate eigenvalues of 7.977 and 2.023 for the sequences A, 
of matrices which we have formed. In particular, these are approximate eigenvalues 
of the matrix Ag which we started with. In fact, the eigenvalues of A, are 8 and 2, 
and clearly, had we continued with the process described above, we would have 
got better approximations to these values. 


The method described above can be expressed as an iterative procedure. Starting 
with a square matrix Ay we form a sequence of matrices A,, r = 0,1,2,... At each 
stage we find the LU decomposition 


A, = L,U, 
and use this to form the next matrix 
AL a U,L,. 


This matrix has the same eigenvalues as A,. 


All calculations here and in 
the rest of this unit are 
worked to the full accuracy of 
the computer. For 
convenience, the numbers 
have been truncated in the 
text. 
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It can be shown that under certain conditions, as r increases, A, approaches an 
upper triangular matrix, with its eigenvalues (and hence those of Ao) in descending 
order of modulus on the main diagonal. The conditions for which the process 
converges are very complicated, and the proof is beyond the scope of this unit. For 
instance, you know that in some cases an LU decomposition does not exist, and 
certainly in those cases we cannot carry out the procedure. 


My advice to you is to first try the procedure and see if it works. If not, try the 
procedure on a matrix of the form A + qi, as you know from Subsection 2.1 that 
this will produce eigenvalues 4, + q, where /; are the eigenvalues you want. 


This method is usually known as the LR method or algorithm (rather than the LU 
method) for this is the way it was referred to when first introduced by Rutishauser 
in 1958. In his original paper, he called the upper triangular matrix R, rather than 
the U adopted in this text. | 


Procedure 3.3: The LR method 


Given a matrix Aj, form a sequence of matrices A, (7 = 0, eae 
using the following procedure: 


1. Find the decomposition 
A, = LU, 


where L, is a lower triangular matrix with ones on its main 
diagonal, and U, is an upper triangular matrix. 


2. Form the new matrix 
A,+1 sas U,L,. 
Do 1-and 2 for r = 0; 12... 


Then, for suitable matrices Ay, and large enough r, A,+, 
converges to an upper triangular matrix, with the eigenvalues in 
order of descending modulus on the diagonal. 


Comment 


If the LR method does not work on Ao, do the same process on 
Ao + ql for some suitable choice of g. In this event, the 
eigenvalues you find will be 4; + q where 4; are the eigenvalues of Ao. 


Exercise 6 


Use three iterations of the method described above to find an approximation for the 
eigenvalues of the matrix : 


6 5 
Ay = 
[Solution on p. 49] 


It should be stressed that this method is for use on a computer and is not meant 
for hand calculation. It is for this reason that there are not many exercises in this 
subsection. You will be given the opportunity to see the LR method in practice 
when you use the computer package EIGSOL, which is described in Section 4. 


This section is only an introduction to a number of methods, called decomposition 
methods, which are currently used for finding eigenvalues. All these methods 
decompose a matrix into two matrices P and Q, and form a new matrix QP. The 
particular example of an LU decomposition was chosen, as it is one of the easiest 
to describe. In the crude form described in this unit, it suffers from slow 
convergence and hence from accumulation of rounding error. The way in which 
the method is amended to speed up convergence is beyond the scope of this unit. 
If you want more information on the numerical methods described in these 
sections, consult The Algebraic Eigenvalue Problem by J. H. Wilkinson (Oxford 
University Press, 1965) which is probably the most comprehensive book to date 
on this topic. 
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Summary of Section 3 
1. If two matrices can be found such that 
A = PQ All the matrices must be 


square. 
we say that A can be decomposed, and the whole process is referred to as a ‘ 


decomposition. 
2. The LU decomposition of a matrix A is 
A = LU 


where L is a lower triangular matrix with ones down the main diagonal, and U an 
upper triangular matrix. Not all square matrices have LU decompositions. 


3. Some important matrix results are: 
(i) AB and BA have the same eigenvalues. 
(ii) A and P~‘AP have the same eigenvalues. 


(iii) The eigenvalues of a diagonal matrix are the elements on the main 
diagonal. 


(iv) The eigenvalues of an upper triangular matrix are the elements on the 
main diagonal. 


(v) Let A be a square matrix, with eigenvalues 1,,12,...,An, and P the matrix 
whose columns are the corresponding eigenvectors x;,X2,...,X, Of A. Then 


Ay 


CS 0 
0 


4. The LR method is a method for finding all the eigenvalues of certain matrices. 
The method is described in Procedure 3.4. It is an example of a decomposition 
method. 


AP = 


ey 


End of section exercise 


Exercise 7 
Do three iterations of the LR algorithm on the matrix 


co 7 5 
ee ees 
Hence make an approximation of the eigenvalues of Ao. 


[Solution on p. 49 | 
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4 The computer package EIGSOL 


4.0 Introduction 
This section describes the computer package EIGSOL. 


EIGSOL is designed to find eigenvalues and eigenvectors of an n x n matrix. 
Again, all the remarks about computer packages in Section 6 of Unit 1 are still 
valid. 


As with other computer packages described in this course, the options are split 
into five main categories. However, before giving a full list of options, we shall 
look at a worked example. 


4.1 A worked example 


In the television programme we used the computer package EIGSOL to find the 
eigenvalue of largest modulus, and the corresponding eigenvector, of the matrix 


6.7 3 2A 
45-31 054, 
23 —5.6 oS 


To do this, only a subset of the options described in Subsection 4.3 is required. 
The following computer dialogue illustrates the options needed for this particular 
problem. The information you have to type in is underlined. The rest of the 
information is supplied by the computer. 


When you have logged on, and obtained the computer package EIGSOL, the 
computer will ask you to name an option. You could then proceed as follows: 


OPTION? 10 
TYPE IN THE SIZE OF THE MATRIX:N=? 3 
OPTION? 11 


TYPE IN THE MATRIX A, ONE ROW AT A TIME, SEPARATING THE 
ELEMENTS BY COMMAS. 


ROW 1? 6.7, 3, 2.1 
ROW 2? 4.5, —5.1, 0.5 
ROW 3? 2.3, —5.6, 9.5 
OPTION? 20 

DIRECT ITERATION 
OPTION? 30 


TYPE IN THE ELEMENTS OF YO, SEPARATING THE ELEMENTS BY 
COMMAS: 


ELEMENTS? 1, 1,1 
OPTION? 33 
SIGNIFICANT FIGURES? 6 
OPTION? 40 

OUTLINE PRINTOUT 


AFTER HOW MANY ITERATIONS DO YOU REQUIRE EACH 
PRINTOUT? 5 


OPTION? SOLVE 


The computer would then respond with the printout given in the television section 
(printout (i) on page 25). 
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4.2 A list of options in EIGSOL 


The options in this package are described under the same general headings as 
those in previous packages. 


Command options 
As before, these options have names, rather than numbers associated with them. 


OPTIONS —_ this option prints the list of options available to you in this 
package. | 
SOLVE —_ this option computes the solution. It is the last option you use 


when solving a problem. Before this is used, all the data and the 
method of solution must have been given to the computer. If any 
of this information is missing, the computer will print an error 
message. 


HELP — help may be obtained if you need it, at any stage that a question 
mark appears. If you type HELP after the computer has printed 
OPTION?, you will be asked questions to help you to decide what 
option to use next. These questions require an answer YES (Y) or 
NO (N). 


If you type HELP while using an option, you will be given help as 
to how to use the option properly. 


LIST — this option prints out the current information held in the 
computer about the current problem. It is advisable to ask for 
LIST before you ask for the problem to be solved. In this way, 
you can check that the computer is holding the correct 
information. 


STOP — this option is your means of exit from the package. You should 
only use this option if 


(i) you want to use another package, or 
(ii) you want to log-off the computer. 


The rest of the options in this package are obtained by typing a number when the 
computer types OPTION? 


Problem options 


10 Specify the size of the matrix 


You will use option 10 to type in the size of the matrix. As all the matrices are 
square in this unit, it is only necessary to type in one number. For example, 
suppose you want to type in a 2 x 2 matrix. When the message 


TYPE IN THE SIZE OF THE MATRIX: N = ? 


appears, just type 2. The maximum size matrix you can use in this package is 
10 x 10. 


11 Enter matrix A 


You will use option 11 to put the matrix A into the computer after you have 
specified the size of the matrix by using option 10. This is done in exactly the 
same way as in the computer package SIMLIN in Unit 9. When the message 


ROW 1 ? 


appears, type in the first row of the matrix, separating the elements by commas. If 
you type too few elements, an error message and an instruction to retype the 
whole row will be printed. If you type too many elements, a warning only will be 
printed. 


12 Edit A 


As in SIMLIN, this option is used to change individual elements of A. When you 
use this option, only the elements you want to alter will be changed. 
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To remind you how to use this option, here is a computer dialogue which will 
change the element a, to be 8. 


OPTION? 12 
CHANGE ELEMENT IN ROW? 2 

AND COLUMN? 3 

TO HAVE THE VALUE? 8 

DO YOU WANT TO CHANGE ANOTHER ELEMENT IN A? NO 


OPTION? 


Had you answered YES to the last question, you would have been asked to specify 
the next element you wanted to change. 


13 Choose a standard problem 


To save you having to type in some of the larger problems at the end of this 
section and in your tutor-marked assignment, there is a set of problems already 
stored in the computer. As in SIMLIN, each stored problem has a name. The 
following computer dialogue will access Exercise 5 (called PROBI) at the end of 
this section. 


OPTION? 13 
PROBLEM NAME? PROB 1 


OPTION? LIST 


It is advisable to print out the standard problem (using LIST) to check that the 
data is what you expected. 


Method options 


20 Direct iteration 


If you use this option, the eigenvalue of largest modulus of A and the 
corresponding eigenvector of A, will be found. To do this, direct iteration as 
described in Procedure 2.4 is used. 


To implement direct iteration you should specify the following parameters: 
(i) the initial vector yo (using Option 30); 


(il) the accuracy which you require in your answer (using Option 33). 


21 Inverse iteration 


If you use this option, the eigenvalue nearest to a specified value p, and the 
corresponding eigenvector of A, will be found. To do this, inverse iteration as 
described in Procedure 2.5 is used. When using this option, you should specify the 
following parameters: 


(i) the initial vector yy (using Option 30); 
(ii) the accuracy which you require in your answer (using Option 33); 
(ii) the value of p (using Option 32). 


22 LR method 


This option should be used when all the eigenvalues of a matrix A + gI are 
required. The method used is the LR method described in Procedure 3.3. To use 
this option, the following values are specified: 


(i) the value of q (using Option 34): 
(ii) a condition for stopping the iteration (using Option 35). 


23 ‘Black Box’ method 


This is an option provided for your convenience. Used in conjunction with 
Options 10 and 11, or with a standard problem (Option 13), this option will 
provide you with the eigenvalues and eigenvectors of any given matrix A, (both 


See the notes on the use of 
these options in the next 
subsection. 
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real and complex). It is called a ‘black box’, simply because you do not know what 
is in it. The method used is a decomposition method, although not the one 
described in this unit. 


You might like to use this option: 


(i) to check the answers you get using the other methods (it replaces the option 
ANSWER in the other packages), or 


(ii) to find quickly any eigenvalues and eigenvectors needed for the problems in 
Units 22 and 24. 


Parameter options 


30 Enter yo 

You will use this option to specify a starting value for yo, both for the direct 
iteration method (Option 20) and the inverse iteration method (Option 21). When 
the message 


ELEMENTS ? 


is printed, type the elements of yo, separating the elements by commas. Typing too 
few elements will produce an error message, and an instruction to retype the whole 
vector. Typing too many elements will just be followed by a warning message. 


If you fail to specify yo, then the first time either Option 20 or 21 is used, yo will 
automatically be set to [1 1 a 1]7. After this, it will retain the value it 
had last time you used it. For instance, if you had set yo to be [1 1 0]’ in 
the previous problem, it will continue to hold this value in the present problem 
unless you change it. 


31 Edit yo 


You will use this option if you want to change individual elements in yo. Suppose 
you wished to change the second element of yo to be 8.4. The computer dialogue 
to do this would go as follows: 


OPTION? 31 

CHANGE ELEMENT IN ROW? 2 

TO HAVE THE VALUE? 8.4 

DO YOU WANT TO CHANGE ANOTHER ELEMENT? NO 
OPTION? 


32 Enter p 


This option is used to specify the value of p in (A — pl)~! in the inverse iteration 
method (Option 21). When the message 


P=? 


is printed, type in the value you choose for p. If you fail to specify p, then the first 
time Option 21 is used, p will be set to zero. After this, p will retain the last value 
you gave it—so be careful. 


33 Specify the accuracy required 


When you use this option in conjunction with Option 20 or 21, it specifies the 
degree of accuracy that you would like y, to have. Suppose you would like two 
successive approximations to agree to five significant figures. Then, when the 
message 

SIGNIFICANT FIGURES ? 


is printed, you type 5. In this case, the iterations will continue until the 
corresponding elements in y,+, and y, agree to five significant figures. 


If you fail to specify this option initially, the iteration will continue until the 
corresponding elements of y,,, and y, agree to three significant figures. After this, 
the accuracy will retain the last value you gave it. You may specify a maximum of 
eight significant figures. 
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34 Enter g 

This option is used to specify the value of q when finding the eigenvalues of 

A + qi using the LR method. Again, if q has never been specified, initially it will 
have the value zero. After this, q will retain the last value you gave it. 


35 Stopping criterion for the LR method 

When you use this option in conjunction with Option 22, it will stop the 
iterations when A, has got sufficiently near upper triangular form. When the 
message 


SUBDIAGONAL ELEMENTS LESS THAN ? 


appears, if you type 0.0001, then the iterations will stop when all the elements 
underneath the main diagonal are less than 0.0001 in magnitude. 


This does not in fact give you much indication as to how accurate your answer 
will be, but obviously, the smaller the value you give this parameter, the greater 
the accuracy you will achieve in your answer. 


If you fail to specify this option initially, then the iteration will continue until all 
the elements underneath the main diagonal are less than 0.0001 in magnitude. 
The smallest value you are allowed to specify in this option is 10~°. 


Print options 


40 Outline printout 
This option gives a printout after each rth iteration. When the message 


AFTER HOW MANY ITERATIONS DO YOU REQUIRE EACH 
PRINTOUT ? 


appears, type the value of r. If you type 5, for example, every fifth iteration in the 
method of your choice will be printed. 


41 Print solution only 


This option will give the solution, and the number of iterations needed to satisfy 
the accuracy restrictions specified; no intermediate results will be printed. If no 
print option is specified this option will automatically be used. 


4.3 Notes on the use of these options 


Use of LIST 


As you have a large choice of method options in this package, you are strongly 
recommended to request a listing of the problem you are about to solve, before 
you ask for the problem to be solved. This is to see that values for yo, p, the 
number of significant figures etc., are all stored correctly before you attempt to 
solve the problem. 


Use of Options 20, 21, 22 


As the number of iterations may be quite large to achieve a sufficiently accurate 
answer, it is advisable to run the option once with ‘solution only’ (Option 41) to 
start with. This will tell you how many iterations were needed to achieve this 
answer, and will help you decide how much detail you would like to see using the 
more detailed printout (Option 40). This advice is given, as printout is very time 
consuming. 


Use of Option 21 

If inverse iteration is used with p = 0, then the routine simply becomes that for 
finding the eigenvalue of smallest modulus of A. 

Use of Option 22 


(i) If all the eigenvalues of A are required, start with gq = 0. Only in the event of 
the LR method not working will you want to try different values of q. 
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(ii) It is particularly important to use ‘solution only’ (Option 41) to begin with, 
for then a suitable message will be printed together with advice to alter the 
value of q if the method is not working. 


(111) If you require the eigenvectors as well as all the eigenvalues, find all the 
eigenvalues using the LR method and then use each eigenvalue in turn as 
your choice for p in the inverse iteration method (Option 21). 


Use of option 23 


This option can only be used with the ‘solution only’ option (Option 41), as the 
intermediate steps are not available for printing. 


4.4 Computer exercises 


You are expected to use the computer package EIGSOL to help you to solve the 
following problems. Some comments on the solutions to Exercises 1 and 3 are 
given on p. 50. 


Exercise 1 


(1) Use direct iteration to obtain the eigenvalue of largest modulus and the corresponding 
eigenvector for the matrix used on television: 


ee 


Use the initial vectors 


3 —2 
(a) o=|2| (b) o=| a 


Can you explain what has happened? 


(ii) Now use direct iteration again, with the initial vectors 


= =, 3 
= oe ee ae (C) Yo = ae 


— 2.00001 
(d) you 1 | 


Explain the results you now obtain in light of the original two results. 


Exercise 2 
Use the computer package to find all the eigenvalues, and corresponding eigenvectors, of 


(i) 2 3 Pao (ii) |61 —2.3 a2. 
10 3 4 1 —3.4 9.7 
3 6 1 See i 11 
In all cases start with an initial vector yo = [1 1 1]7, and continue until two 


consecutive iterates agree to six significant figures. Proceed as follows: 
(a) Use direct iteration to obtain the eigenvalue of largest modulus. 
(b) Use inverse iteration to find the eigenvalue of smallest modulus. 


(c) Suggest a suitable value of p, for use in the inverse iteration method to find the 
remaining eigenvalue, and hence find the remaining eigenvalue. 


Check your results using Option 23: 


Exercise 3 


(1) Use the LR method to find the eigenvalues of the matrices in Exercise 2. Specify six 
significant figure accuracy using Option 35. 


Your answers may differ slightly from those in Exercise 2. Can you explain this? 
(11) Use the eigenvalues you have just found as suitable values of p in the inverse iteration 
method to obtain the corresponding eigenvectors. , 


Exercise 4 
(1) Check that the LR method fails for the matrix 
0 Z 3 
A=i2 3 0 |. 
1 7 10 
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(ii) Use the LR method to find the eigenvalues of a matrix A + ql, with a suitable choice 
of q. Hence find the eigenvalues of A. 


Exercise 5 


(To save you time, the coefficients of this problem are already stored in the computer, under 
the name PROBI. Use Option 13 to obtain PROBI.) 


Find the eigenvalue 
(i) of largest modulus 
(11) of least modulus 
of the matrix 
12.4 if a ie ae 
17 10.8 —-1.1 oH | 3.1 
A=| 35 —1.1 7.8 —3.0 0.1 |. 
1.1 2.7 —3.0 8.4 1.4 
—2.7 31 0.1 1.4 9.6 
(iii) Use trial values of p in the inverse iteration method to find the other eigenvalues of A. 


(iv) Use the LR method to check the eigenvalues which you have found. 


Summary of Section 4 


This section describes the options in the computer package EIGSOL. This package 
implements the methods described in Sections 2 and 3, all of which find 
eigenvalues of a matrix A. 


The following lists of options are the minimum needed to implement these 
methods. 


Description of option LR method 


Direct iteration Inverse iteration 


Specify the size of A 10 10 10 
Enter A 11 1 11 
Method to be used 20 21 22 
Specify p not applicable x not applicable 
Specify g not applicable not applicable 34 
Print solution only 41 41 41 
Advisable to list 

before solving Last LIST LIST 
Solve the problem SOLVE SOLVE SOLVE 


5 End of unit test 


The first four exercises cover routine material from the unit and you should make 
sure you can tackle them with confidence. The remaining exercises may require a 
little more thought. 


The ‘first three exercises all refer to the matrix 


: 7 14 
i 


Exercise 1 
Find, directly, the eigenvalues and eigenvectors of A. 
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Exercise 2 


Use suitable iterative methods to find approximations to the following selected eigenvalues 
and corresponding eigenvectors of A: 


(i) the eigenvalue of largest modulus, 

(ii) the eigenvalue of smallest modulus, 

(iii) the eigenvalue closest to 12. 

In all cases, use a starting vector yo = [1 1]7, and continue the iterations until the 


elements in two iterates y,,, and y, agree to two significant figures. 


Exercise 3 
Use three iterations of the LR method to make an approximation of the eigenvalues of A. 


Exercise 4 


Given that A is an n x n matrix with distinct eigenvalues 4,,/,,...,1,, and corresponding 
eigenvectors X;,X>,...,X,, find the eigenvalues and eigenvectors of the following matrices in 
terms of the eigenvalues and eigenvectors of A: 


(i) A° 

(ii) (A+ pl)~*(A — pl) 

Exercise 5 

In certain cases, symmetric matrices have a decomposition LL’, where L is a lower 


1a. = = 
triangular matrix. For instance, the 2 x 2 matrix b has the decomposition 
c 


Ja O 1f a b 


es 
Zz [actA |} 0 deta] 
Ja a a 


so long as a > 0 and det A > 0. This type of decomposition is called the Cholesky 
decomposition. 


(i) Show that L’L is also a symmetric matrix. 


(ii) Find the Cholesky decomposition of the matrix 


Exercise 6 
Consider the following iterative scheme: 


Given a symmetric matrix Ao 

1. Find the Cholesky decomposition 
A= LI. 

2. Form the new symmetric matrix 
An =i 

Do 1 and.2 for? = 0,1,2... 


In this way a series of symmetric matrices Ap, A;, Az,... are formed. 
(i) Show that the sequence of matrices formed all have the same eigenvalues. 


(ii) Find A,, A, and Az; for 


cS 9 8 
"1s. 3) 
Hence, find an approximation for the eigenvalues of Ao. 


Exercise 7 


Suppose that a given matrix A has both an LU decomposition and a Cholesky 
decomposition, and that both the LR method and the method described in Exercise 6 work. 
Which method would you use to find an approximation to the eigenvalues of A, and why? 


[Solutions to Exercises 1 to 7 on pp. 50-52 | 


There are, of course, only two 
eigenvalues and so the 
iterations in (iii) will converge 
to one of the eigenvalues 
found in (i) and (ii). The 
purpose of (iii) is to give you 
extra practice. 
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Appendix: Solutions to the exercises 


Solutions to the exercises in Section 1 


1. (i) For the matrix 


ee & , 
3 ‘ 
the condition det (A — AI) = 0 gives 
T—A 3 
3 T—A 
Hence we obtain the equation 
(7-—Aay’r—3% =0 
or (4—A)(10 —A)=0 


and so the eigenvalues of A are 4 and 10. 


[= 


(ii) For the matrix 
3 2 
As 
hod 
det(A — AI) = 0 gives 


_ 2 
P A — 


1 4d 
Hence we obtain the equation 
(33 —A)(4—A4)-2=0. 
This gives us 
A*—74+10 =0 
or (A—5)(A —2)=0. 
Hence, the eigenvalues of A are 5 and 2. 
2. (i) The equations determining the eigenvectors are 
(7 — A)x, + 3x, = 0. 
3x, + (7 — A)x, = 0. 
We shall consider each eigenvalue in turn. 


Case 1 = 4: in this case the equations become 


3x, + 3x, = 0. (twice) 
Thus we get a solution x; = —k, x, =k. So, an eigenvector 
corresponding to 4 = 4 is - 
Case 4 = 10: in this case the equations become 
— 3x, + 3x, =0 
3x, — 3x, = 0. 


(Note that these equations provide us with the same 
information. ) 


These have a solution of the form x, = k, x, = k. So, an 


; ; ; i 
eigenvector corresponding to A = 10 is it 


(ii) The equations determining the eigenvectors are 
(3 — A)x, + 2x, =0 
xX, + (4-—A)x, = 0. 
We shall consider each eigenvalue in turn. 


Case A = 2: in this case the equations become 


xX, + 2x, =0 
Xj + 2X = 0. 
Putting x, = k, we get x, = —%4k. So an eigenvector 


corresponding toA=2is| _ 4 |: 


Case 4 = 5: in this case the equations become 
—2x, + 2x, =0 
ty = X> = 0: 
These equations have a solution of the form x, = k, x; =k. 


1 
So an eigenvector corresponding to 1 = 5 is i 
3. The condition det(A — AI) = 0 gives us 
1-—A ae 6 
0 4—} 4 j=. 
0 0 2—A 


Hence (1 — 4)(4 — 4)(2 — A) = 0, and the eigenvalues are 1, 
4, and 2. To find the corresponding eigenvectors, we solve 
the equation 


1-—A 2 6 x4 
0 4—j) 4 X2 =: &, 
0 0 2-—A/ X3 


putting 4 equal to each eigenvalue in turn. 


Case A = 1: in this case we obtain the equations 


2X2 + 6x3 = 4) 
3x5 + 4x, 4) 
xy =. 


This gives us x3 = 0, x, = 0. As x, is not involved in these 
equations, it can take any value k. 


Hence [1 0 
ee 


Case 4 = 4: in this case we obtain the equations 
— 3x, > 2X5 + 6x3 =) 
+3. >= 0 


0]" is an eigenvector corresponding to 


(twice). 
Thus 

— 3x, + 2x, = 0, 
and putting x, =k, we get x, = 3k. 


So [1 3 0]? (or [2 3 
corresponding to 4 = 4. 


0]") is an eigenvector 


Case 4 = 2: in this case we obtain the equations 
—X4 ee 2X2 + 6x3 a 
2x, + 4x3 = 0. 


So putting x, = k we find that x, = —2k and x, = 2k. So 
[2 —2 1]" is an eigenvector corresponding to 4 = 2. 


4. From Example 7 we know that the eigenvalues of B are 
1 + 2i and 1 — 2i. To find the corresponding eigenvectors, we 
solve (B — AI)x = 0, that is 


os —1 |e \=0 
4 1—A]| x, : 
putting A equal to each eigenvalue in turn. 
Case A = 1 + 2i: in this case we obtain the equations 
—2ix; — x,=0 
4x, — 2ix, = 0. 


(Note that the second equation is simply 2i times the first 
one.) From either equation, we obtain solutions x, = ik/2, 
X27 = k. 
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i/2 
Hence, an eigenvector corresponding to 1 + 2i is ¥ | 


Case 4 = 1 — 2i: in this case we obtain the equations 
2ix,; — x, =0 
4x, + 2ix, = 0. 
From either equation, we obtain solutions of the form 
X, = —ik/2, x2 = k. 


Bans 
Hence, an eigenvector corresponding to 1 — 2i is : | 


5. (i) The characteristic equation is given by 


1-—A 2 

ee 
that is 

(a —1)4=—1)43 6 
or A* —54+6=0 


Le. (A — 2)(A — 3) = 0. 


So the eigenvalues are 2 and 3. To obtain the eigenvectors, 
we solve the equation 


1-—A 2 ¥ 
oo 
putting A equal to each eigenvalue in turn. 
Case 4 = 2: we obtain the equations 
—x, + 2x, =0 
—x, + 2x, =0 


which have solutions of the form x, = 2k, x2 = k. Hence 
[2 1]’ is an eigenvector corresponding to the 
eigenvalue 2. 


Case 1 = 3: we obtain the equations 
—2x; Se 2X> == () 
—X,;+ x,=0 


which have solutions of the form x, = k, x. = k. Hence 
[1 1]" is an eigenvector corresponding to the 
eigenvalue 3. 


(ii) The characteristic equation is given by 


4—; S * 
5 TA 
that is 
(4 —A\(7 —4)— 10 = 0 
or A? — 114 + 18 =0 
i.e. (A —9)(A — 2) = 0. 


So, the eigenvalues are 9 and 2. To obtain the eigenvectors, 
we solve the equation 


4— A e X41 = 0 
5 Vicdtt ae LS 
putting 2 equal to each eigenvalue in turn. 


Case 4 = 9: we obtain the equations 


—5x, + 2x, =0 
3X4 = 2X2 = 0 
which have solutions of the form x, = k, x; = 3k. Hence, 
[2 5]' is an eigenvector corresponding to the 


eigenvalue 9. 

Case 4 = 2: we obtain the equations 
2x, + 2x, =0 
5x, + 5x, =0 
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which have solutions of the form x, = k, x; = —k. Hence 

[1 —1]’ is an eigenvector corresponding to the 

eigenvalue 2. 

(iii) The characteristic equation is given by 

1-A —1 
1 1-A 


: 


that is, 
(1—A)*>+1=0, 
or (1—-A+i)1-—A-i) =0. 


so, the eigenvalues are 1 + i and 1 — i. To obtain the 
eigenvectors, we solve the equations 


is =i 41x 

alhe |= 
putting A equal to each eigenvalue in turn. 
Case 4 = 1 + i: we obtain the equations 

—ix,; — xX, =90 

x, — ix, = 90, 

which have solutions of the form x, = ik, x, = k. Hence, 
[i 1]’ is an eigenvector corresponding to the eigenvalue 
L+t 
Case 4 = 1 — i: we obtain the equations 

ix; — X2 =0 

x, + ix, =0, 
which have solutions of the form x, = —ik, x, =k. Hence 


[-—i 1]’ is an eigenvector corresponding to the 
eigenvalue | — i. 


(iv) The characteristic equation is given by 


1-—A 0 —1 
1 2—A 1 j= 0. 
2 2 3-—A 
Expanding the determinant, we obtain the equation 
2—A 1 1 2-—A 
a-a)"5 sal-/ pS 
This gives 


(1 — a)(2 — A)(3 — A) -2(1 — A) — 2-22 —4)) =0 


which simplifies to give 
(1 — A)(2 — A)(3 — A) = 0. 


So the eigenvalues are 1, 2 and 3. To obtain the eigenvectors, 
we solve the equation 


1-—A 0 —| X4 
1 2-A 1 X2 = 0 
2 2 3-—A X3. 


putting A equal to each eigenvalue in turn. 
Case 4 = 1: we obtain the equations 
—x,=0 
Xi + X2+ x3 =0 
2x, + 2x2 + 2x3 =0 


which have solutions of the form x; = —k, x2 =k, x3 =0. 
Hence [—1 1 0]’ is an eigenvector corresponding to 
the eigenvalue 1. 


Case 4 = 2: we obtain the equations 
—X, —x,=0 
x; +x3=0 
2x, + 2x, + x3 =0 
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which have solutions of the form x3, = k, x, = —k, and 
x2 = k/2. 


Hence, [ —2 1 
the eigenvalue 2. 


2]’ is an eigenvector corresponding to 


Case 4 = 3: we obtain the equations 


—2x, —x;=0 
Xi — X2+x3=0 
2x; + 2x, = ( 
which have solutions of the form x3 = k, x, = —k/2, 
Xo KZ. 


Hence [—1 1 
the eigenvalue 3. 


2]' is an eigenvector corresponding to 


(v) The characteristic equation is given by 


3-4 z 2 
2 2-A QO =f 
2 0 4—} 
Expanding the determinant, we obtain the equation 
2-—A 0 2 0 
@ 2) 0 a 2 hep 
2 2-—A 
of han 
Hence 
(3 — A)(2 — A)(4 — A) — 4(44 — A) — 422 — A) =0 
sO (3 — A)(2 — A)(4 — A) — 833 — A) =0 
ie. (3 — A)(A* — 64) = 0. 


Thus, we obtain the eigenvalues 0, 3, and 6. To find the 
corresponding eigenvectors we solve 


3-—A Z 2 x; 
2 2-—A 0 Xx, 1=0 
Z 0 4-j X3 


putting 2 equal to each eigenvalue in turn. 
Case 4 = 0: we obtain the equations 


3x, + 2x, + 2x3 =0 


ox “f- 2x> == Q 
2X5 + 4x, = 0, 
which have solutions of the form x3 = k, x; = —2k, x, = 2k. 


Hence, [—2 P 
the eigenvalue 0. 


1]’ is an eigenvector corresponding to 


Case 4 = 3: we obtain the equations 


2x, + 2x3; =0 
2x; —-— X2 ==) 
oes + x,=0 
which have solutions of the form x3 = k, x, = —k, 
x, = —k/2. 
Hence, [1 2  -—2]' is an eigenvector corresponding to 


the eigenvalue 3. 


Case 4 = 6: we obtain the equations 


—3x, + 2x2 + 2x3 <= () 


2x = 4x, = 
2X1 — 2x; =0 
which have solutions of the form x, = k, x, =k, x, = k/2. 


Hence, [2 1 
the eigenvalue 6. 


2]’ is an eigenvector corresponding to 
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Solutions to the exercises in Section 2 


1. We have seen that for any number gq, the matrix (A + ql) 


has eigenvalues 4, + q, A.+4q, ..-, A, =, SO 
(i) A-+5I has eigenvalues 1,+ 5, 4,+5, ..., A, +5. 
(ii) A — 2I has eigenvalues 14, -—2, A,—2, ..., A, —2. 
(iii) A — gl has eigenvalues 1, —q, 4,—q, ..., A, —@. 
2. We know 

Ax; = 4;X; tia 122 at 


Subtracting px; (or equivalently pIx;) from both sides of this 
equation, we get 

Ax; — pIx; = 1;X; — px; for i= 1.2.13. 
so (A — pl)x; = (4: — p)X: 


Since (A — pI) is non-singular, we can multiply both sides of 
this equation by (A — pI)~/ to obtain 


x; = (A; — p)(A— pl)” *x; 


for i= 162... En. 


fat te 1 2. 


or (A -— pl)~'x; = X; for i= 1,2,...,n. 
A; — p 
Thus (A — pI)~* has the set of eigenvalues 
1 1 1 

So, *' 7, and has the same cigenvectors 

oP eee Ap 
as A. 
3. Step 1: 


oe ee 
w-[s alll-La} 
Dividing by 12 gives y, = [0.5 iy. 
Step 2: 


1-2 4 

Rebs als |*lest 
Dividing by 9.5 gives y, = [0.421 1)". 
Step 3: 


,_[4  2]f0.421]__ [3.684 
Pe ito 7b hea 005) 
Dividing by 9.105 gives y, = [0.405 iY. 


Thus we obtain an estimate of 9.105 for the eigenvalue of 
largest modulus and a corresponding eigenvector 
[0.405 Lr. 


These compare well with the actual value 9 and [0.4 =. 


7 3 eee | 5 -3 
4. Given A =| { sf then A =a _g - 


Using the scheme y,., = A” 'y, with scaling, and an initial 
vector yg = [1 1]’ we obtain 


Step 1: 
rs -3y 0.182 
ase alos Aa 2 faa | 
Dividing by 0.182 we obtain y, = [1 —0.5]". 
Step 2: 


— a 6 a 0.591 
—— al os sere) = eA 
Dividing by —1.045 we obtain y, = [-0.566 1]. 
Step 3: 


240 & —330-oset (f —05%0 
oT ok Td ee tee 


Dividing by 1.047 we obtain y; = [ —0.506 ee 
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Thus, the estimated eigenvalue of smallest modulus for A is 
1/1.047, i.e. 0.955, with a corresponding estimated eigenvector 
[ — 0.506 1}’. Although this is not a bad approximation 
to the actual eigenvalue 1 and eigenvector [—0.5 LY, 
until the values of y, settle down a bit more it would not be 
very sensible to stop at this point in practice. 


1 0 1 
| and p= 2, then A ~ pt =| 


z 
5. Given A=| 4 1 5 


and so 


ee 
(A — pl) “ls wd 


Using the scheme y,,, = (A — pI)” *y, with scaling, and an 
initial vector yo = [1 0]? we obtain 


Step 1: 


cor Pag 28 gy 1T oP 
e=[T ollol-L a 
Dividing by —2 we obtain y, = [1 —0.5]’. 
Step 2: 


= to 25 
mL t ollosl-[ 7] 
Dividing by —2.5 we obtain y, = [1 —0.4]’. 

Step 3: 


1-2 te 
4 4 0) 
Dividing by —2.4 we obtain y; = [1 —0417]’. 


Hence, an estimate for the eigenvalue of A nearest to 2 is 


1 
oy + 2 or 1.583, and an approximation to the 


corresponding eigenvector is [1 —0.417]". The actual 
eigenvalues are 4.412 and 1.588 (to three decimal places). 
Although our approximation was not too bad, it would have 
been advisable to continue, as the process had not settled 
down sufficiently to quote this result with confidence. 


6. (i) The direct iteration scheme y,,, = Ay,, with scaling 
gives 
Step 1: 
330 40] =p — 39 
ee ose 4 e 24: 
Dividing by —39 givesy, = [1 0.513)’. 
Step 2: 
ee ee ee — 18.487 
= a lane Z yey | 
Dividing by — 18.487 gives y, = [1 0.499’. 
Step 3: 


_ [39 40]f 1 ]_[-19.028 
Y3=| 959 21 |10499| | —9.515 | 
Dividing by — 19.028 gives y; = [1 0.500)’. 


Hence an estimate for the eigenvalue of largest modulus is 
— 19.028 and an estimate for a corresponding eigenvector is 
(1 0.5)". 


Pee —~ 
(ii) With p=2,—pl=| ~ 9 a ees 
Ee a 
eat =a 30 —4i} 


Using the inverse iteration scheme y,+; = (A — pl)~*y,, with 
scaling gives 
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Step 1: 
1719 -401[1] [0.905 
tie lo ee = bese 
Dividing by 0.952 gives y, = [0.95 tT: 
Step 2: 
| 1719 4017095] _ [—1.045 
ae sil 20 we | 1 . : ied 
Dividing by — 1.048 gives y, = [0.998 ee 
Step 3: 


pa 0.998 | | —1.0021 
— sil 2 sy 1 |- Be 
Dividing by — 1.0022 gives y3 = [1.000 i. 
Hence an estimate for the eigenvalue of A nearest to 2 is 


i 


— 1.0022 
corresponding eigenvector. 


+ 2 = 1.002, and [1 1]? is an estimate for the 


Solutions to the exercises in Section 3 


1. We require 


fale ile 2 


Hence 
b= 6, e= 5, 
ab = 4, ac+d=5. 
So a = 4, and d = 3. 


Thus we obtain the decomposition 


Esl ile 3] 


2. We require 


1 yi 3 1 0 O|| d e 5 
1 6 4\/=J]a 1 0}| 0 g h 
2 1 4 b C 1||.0 0 j 
d e Fa 
=| ad ae+g af +h 
bd be + cg bf + ch+j 


Equating coefficients, we obtain: 

First row: d= 1, 2 = 2, f = 3. 

First column (remainder!): ad = 1, bd = 2, 
which gives a = 1 and b = 2. 

Second row: ae+g=6,af+h=4, 
which gives g = 4 andh = 1. 


Second column: be + cg = 1, 


which gives c = —3. 
Third row: bf + ch + j = 4, 
which gives j = —3. 


Thus the decomposition is 


1 3 3 i 0 Oj} 1 2 3 
1 6 4;=)1 1 0}) 0 = 1}. 
2 1 S 2 -7 1} [0 0 -2 


We can check this by multiplying out the right-hand side. 
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Note that there is more than one order in which you can 
equate the coefficients. The row and column approach used 
here is quite a standard one. You can equally well find the 
coefficients row by row, or column by column. 


3. First consider AB 


e-( g 


The eigenvalues of AB are given by det(AB — AI) = 0. That 
is 


1-—A 4 
3 , {no 
giving 
(i= HG—e Se 
or A? — 34-10 =0 
ie. G30 +2 = 


Thus we obtain eigenvalues 5 and —2 for AB. 


Now consider BA 


ei ct sift @) Pe 4 
VQ 1004 24 a 


The eigenvalues of BA are given by det(BA — AI) = 0. That 
1S 


13-1 4 {4 
—30 -10-A 
giving 
(13 — 4)(—10 — 4) + 120 =0 
or P= 
ice. (A—5)A+2)=0 


and hence, once again we obtain eigenvalues 5 and —2. 

So, in this case AB and BA have the same eigenvalues. 

Note: It would have been sufficient to show that AB and BA 
had the same characteristic equations. 


2 0 
4. (i) Using the result of Theorem 3, r | has 
eigenvalues 2 and 3. 


7 2 1 
(ii) Using the result of Theorem 4,| 0 8 5 | has 
eigenvalues 7, 8 and 9. 0 0 9 


5. The matrix A has eigenvalues given by 
? —i 3 
3 2-—A 
or (2—A)* —37=0. 
This gives eigenvalues 5 and —1. Substituting these values 
into (A — AI)x = 0, we obtain an eigenvector [1 ry 


corresponding to the eigenvalue 5, and [1 —1]? 
corresponding to the eigenvalue — 1. 


=0 


Thus Theorem 5 gives us that 


1 1 5 0 
> eel 8 


and you can check in this case that 
P-'AP =D. 


The matrix P is not unique, for it depends on the choice of 
eigenvectors. You should check that P~1AP = D for your P. 
There are only two possible answers for D—the one above 


Or 
— <= o 
2 


if the columns of P are reversed. 
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In other words, D must be a diagonal matrix with the 
eigenvalues of A on the diagonal. 


6. We shall apply Procedure 3.3 to the matrix 


6 5 
A, = 
le | 
recording all numbers to three decimal places. 


First iteration: 


gies 1 0116 5 
eee i ee eee, 


A, =Usto =[° 5 l 1 é 
0 1.667||0.667 1 
9.333 5 
-|* ‘111 pe 


1 0 1[ 9.333 5 
A, = L,U, = 
ge ae = Alf 0 7 


9.333 5 1 0 
A, =U,L, = 
Sa 0 i. | 


S 9.929 5 
40,198 = 2071 
Third iteration: 


0 | | 9.929 5 
Aa Slabs i. 613 , 0 a 
9.929 5 1 0 
i ee -| 0 ee ae N 
9.993 5 
x he fa 
Hence, an estimate for the eigenvalues of A would be 9.993 
and 1.007 after 3 iterations. 


You can check that the eigenvalues are in fact 10 and 1. 


7. We start with the matrix 
7 5 

roe X 7 

First iteration: 
0 
oon - 429 1 ks v3 aay 
7 
oe E 2. oe Fe 429 4 
143 5 
- 1.224 Set 

Second iteration: 


A, =L,U, = 


Third iteration: 


0 || 9.813 5 
= | 
ue =| 00 0299 ee 
9.813 5 
A; = ULL, = ce | 
| 9962 — 
~ 10.0609 —- 2.038 | 
Hence, an approximation for the eigenvalues of Ag after 
three iterations would be 9 962 and 2.038. 
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Solutions to the exercises in Section 4 


The answers to most of the exercises can be checked using 
the ‘black box’ routine (Option 23). However, Exercises | 
and 3 require the additional explanations given below. 


1. (@ Starting with the initial vector yo = [3 1]? (Case 
(a)) the scheme produces the eigenvalue (and corresponding 
eigenvector) of largest modulus, as expected. However, 
starting with the initial vector yo = [—2 
scheme produces the eigenvalue of smallest modulus. This 
unexpected result can be explained by observing that the 
scheme stipulates that the initial vector yo must not be an 
eigenvector. In Case (b) the initial vector was chosen to be 
Yo=[-2 1]’, which is an eigenvector corresponding to 
the eigenvalue of smallest modulus. 


You might well ask how we can avoid choosing yo to be an 
eigenvector. If the scheme reaches the solution in only one or 
two iterations, you have reason to be very suspicious! 


(ii) The initial vectors yy—however close they are to an 
eigenvector corresponding to the eigenvalue of smallest 
modulus—are not in fact eigenvectors and do, as the 
geometry on the television programme suggests, get dragged 
round to the dominant eigenvector. 


3. (i) The answers differ slightly from those in Exercise 2 

because the LR method stops when the magnitude of all the 
elements below the diagonal in A, are less than the tolerance 
you stated. This does not in fact tell you how accurate your 

answer is going to be. It just ensures reasonable accuracy in 

the answer you give. 


Solutions to the end of unit test 
1. det(A — AI) = 0 gives 


7T-—A 14 
3 5 n8 
Hence, 
J 8-6-0 
or A? — 1544+ 14 =0 


Le. (A — 14)(4 — 1) =0. 


So, the eigenvalues of A are 1 and 14. The corresponding 
eigenvectors can be found by solving (A — AI)x = 0, putting 
A equal to each eigenvalue in turn. 


Case 4 = 1: we get the equations 


6x, + 14x, st 
3x1 + 7X2 <= 0. 
Thus we obtain a solution x, = k, x, = —4k, and hence an 


eigenvector [—7 3]’ corresponding to the eigenvalue 1. 


Case 4 = 14: we get the equations 
— 7x, + 14x, =0 
3X1 net 6x, oe). 


Thus we obtain a solution x, = k, x, = 2k, and hence an 
eigenvector [2 1]? corresponding to the eigenvalue 14. 


2. (i) To find the eigenvalue of largest modulus, we use 
direct iteration with scaling. 


Step 1: 


< to ee 
vey ‘li-ti} = 


So, dividing by 21 we obtain y, = [1 0.524)". 


Step 2: 


7 1417 1 14.333 
ie e = 14.333. 
de t | ee ceo %2 


1]7 (Case (b)) the 
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So, dividing by 14.333 we obtain y, = [1 0.502 ]’. 


Step 3: 


ear ae ee 
9715 gos) | Fone 


So, dividing by 14.023 we obtain y3 = [1 0.500]7. 


As y, and y; agree to 2 significant figures, we shall stop the 
iteration. Hence, to two significant figures the eigenvalue of 
largest modulus is 14, and a corresponding eigenvector is 

[1 0.50)". 


(ii) To find the eigenvalue of smallest modulus, we use 
inverse iteration y,,, = A ‘y,, where 


2 4g 
pie ag 7 


es: a, = —0.429. 
So, dividing by —0.429 we obtain y, = [1 —0.667]’. 
Step 2: 
4 fF Be 
pat al 3 | ceed 
of 28 one 
So, dividing by 1.238 we obtain y, = [1 —0.442]’. 
Step 3: 
1 8 —14 1 
ee al 3 < eae 
_| ae ata 
So, dividing by 1.014 we obtain y; = [1 —0.430]’. 
Step 4: 
* 4) ga 1 
sags al 3 - ae 
1.001 
a | at gf = LOOT. 
So, dividing by 1.001 we obtain y, = [1 —0.430]’. 


Now, y. and y; agree to 2 significant figures and so we shall 
stop the iteration. Hence, an approximation for the 
eigenvalue of smallest modulus of A is 1/1.001 = 1.0 (to 2 
significant figures). The corresponding eigenvector is 

[1 —0.430]’. 


—-5 14 
(iii) a-121=| 3 SO 
1|4 
A — 121)! =— 
( ) | 3 5 


Yr+1 = (A — 121) “*y, gives 


| Inverse iteration 


Step 1: 


oo ok [40 ai 
as. 5a 


0.818 
= a Z,; = 0.818. 


So, dividing by 0.818 we obtain y, = [1 0.444)". 
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Step 2: 


ste ta 
q2 993) 2 Hl 0.444 


0.465 
-| i a“ = 0.465. 


0.237 
So, dividing by 0.465 we obtain y, = [1 — 0.511]. 
Step 3: 
iaest [Anse eet 9 
Joe al 3 ee) 
= bags a3 = 0.507. 
So, dividing by 0.507 we obtain y; = [1 _—_ 0.498]. 
Step 4: 
<< ee Ga 4 
ee al = lle 
(a. a 
So, dividing by 0.499 we obtain y, = [1 0.500". 


Now, y3 and y, agree to 2 significant figures, and so we shall 
stop the iteration. Hence, an approximation for the 


eigenvalue nearest to 12 is + 12 = 14 (to 2 significant 


0.499 
figures). The corresponding eigenvalue is [1 Osi. 


3. First iteration: 


An = £8 " 
-t 2 


7 0 
oe i 7 429 : 
13.000 
-| 0.857 = 
Second iteration: 
13 14 
eo e 0659 a 0 ee 


13 1 0 
0 1.077 | |.0.0659 1 


13.923 14 
0.07100 1.077 


13.923 14 
Fe Fe: 00510 | 0 one 
13.923 1 0 
ce. “5 1 pa ee ‘ 
13.994 14 
= irae fee 
Thus, approximations for the eigenvalues of Ay are 13.994 
and 1.006—or 14 and 1.0 to two significant figures. 


Third iteration: 


4. We know that 
Ax; = 1;X; gt os Se eee (1) 
(i) Multiplying Equation (1) by A, we get 
A’x; = 1,Ax; 
= A,(A;X;) (using (1)). 
So A°x, = A7x; (2) 
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Multiplying Equation (2) by A, we get 


A°x; < 42 Ax; 
= A?(A;x;) (from (1)). 
So A = iF: 


Thus, the eigenvalues of A° are 1}, A3,...,A3, and the 
corresponding eigenvectors are the same as those of A, i.e. 
ore Ceres, eo 


(ii) Adding pIx; (or equivalently px;) to both sides of (1), we 
get 
(A + pI)x; = (A; + p)xi. (3) 


Similarly, subtracting pIx; (or equivalently px;) from both 
sides of (1), we get 


(A — pI)x; = (4; — p)x;. (4) 
| Re-arranging (3) and (4) we get 
a ea (Ae 
A, + p A; — p 
Equating these two, we get 
(A — pI)x; (A + pI)x; 
4p Ate 
Left-multiplying this relation by (A + pI)~‘(A; — p) we get 
= 4: — P 
(A + pl) (A — pl)x; = Ey a 


for 


. A; ae 
So, the eigenvalues of (A + pI)~*(A — pl) are P 
P= 4,2... 0 


The eigenvectors are the same as those of A, ie. X,,X>,...,X,. 
(You can probably see by now the connection between the 
eigenvalues of A and those of various ‘functions’ of A.) 


(i) 
Ja = 0 


° /detA 


which is a symmetric matrix. 


1 
—detA 
a 


(ii) From the given decomposition: 
ke | = 3 0 I aprd 
8 9 2.667 1.374 || 0 1.374 
6. (i) We have 
Ao = LoLé 
Ae = Esko. 


We know from Theorem 1 of Section 3 that, for any square 
matrices A and B, AB and BA have the same eigenvalues 
and so Ay and A, have the same eigenvalues. Similarly, so 
do A, and A, etc. Thus the sequence Ao, A,,A>,... all have 
the same eigenvalues. 


(ii) To cut down the working, we shall use the result of 


b 
Exercise 5(i) which tells us that if A, = e then 
c 


b 
° det 


1 
—detA 
a 


First iteration: for Ao, a = 9, b = 8, detA = 17, and so 


Second iteration: a = 16.111, b = 3.6650, det A = 16.9998, and 


3.6650 x ./16.9998 
~ | 3.6650 x ./16.9998 
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Third iteration: a = 16.9447, b = 0.93794, detA = 16.9998, 
and so 


0.93794)? ./ 16.9998 
16.9447 + ( ) 0.93794 x ——_—__— 


Pare 16.9447 16.9447 
ee aes /16.9998 16.9998 
* "16.9447 16.9447 


16.9966 0.2282 
0.2282 1.0033 | 


As we know that Ay, A;, A> and A; all have the same 
eigenvalues, we can get an approximation from A; for the 
eigenvalues of Ay of 16.9966 and 1.0033—or 17 and 1.0 to 
two significant figures. 


7. Given that both methods work, it would be advisable to 
choose the LL’? decomposition, simply because there are less 
elements to find—and hence, less work involved. 


re 
q : yf 
) + 
ae | fea 
{ i r ( 
q 
j t 
; ‘ 
' ’ \ . 
t | 
4 “& 
} : | | 

| " % q 
aes | ib 

; ‘ : | 

4 : 

4 | | 

guar” 
P : 

ar an as ‘ 
{ .  & wy = fi 7 a : Ag - 


an 


vn 
U 


y i 


